Challenge: What is the future selling price of a home?

A home is often the largest and most expensive purchase a person makes in his or her lifetime. Ensuring homeowners have a trusted way to monitor this asset is incredibly important. In this competition, students are required to develop a full-fledged approach to make predictions about the future sale prices of homes. A full-fledged approach constist, at least, in the following steps:

Descriptive statistics about the data
Data cleaning and pre-processing
Defining a modeling approach to the problem
Build such a statistical model
Validate the outcome of the model

Now, should you ask a home buyer to describe their dream house, they probably wouldn't begin with describing features such as the height of the basement ceiling or the proximity to a railroad. As you will see, the dataset we use in this competition proves that many more features influence price negotiations than the number of bedrooms or a white-picket fence. With 79 explanatory variables describing (almost) every aspect of residential homes in a small city in the US, this competition challenges you to predict the final price of each home.

#this is just to remove the warnings
import warnings
warnings.filterwarnings('ignore')
warnings.filterwarnings('ignore', category=DeprecationWarning)

#installing some libraries
!pip install pandas_profiling
!pip install xgboost

Collecting pandas_profiling
  Downloading https://files.pythonhosted.org/packages/a7/7c/84f15ee705793a3cdd43bc65e6166d65d36f743b815ea517b02582989533/pandas_profiling-1.4.1-py2.py3-none-any.whl
Requirement already satisfied: jinja2>=2.8 in /opt/conda/lib/python3.6/site-packages (from pandas_profiling)
Requirement already satisfied: six>=1.9 in /opt/conda/lib/python3.6/site-packages (from pandas_profiling)
Requirement already satisfied: pandas>=0.19 in /opt/conda/lib/python3.6/site-packages (from pandas_profiling)
Requirement already satisfied: matplotlib>=1.4 in /opt/conda/lib/python3.6/site-packages (from pandas_profiling)
Requirement already satisfied: MarkupSafe>=0.23 in /opt/conda/lib/python3.6/site-packages (from jinja2>=2.8->pandas_profiling)
Requirement already satisfied: python-dateutil>=2 in /opt/conda/lib/python3.6/site-packages (from pandas>=0.19->pandas_profiling)
Requirement already satisfied: pytz>=2011k in /opt/conda/lib/python3.6/site-packages (from pandas>=0.19->pandas_profiling)
Requirement already satisfied: numpy>=1.9.0 in /opt/conda/lib/python3.6/site-packages (from pandas>=0.19->pandas_profiling)
Requirement already satisfied: cycler>=0.10 in /opt/conda/lib/python3.6/site-packages (from matplotlib>=1.4->pandas_profiling)
Requirement already satisfied: pyparsing!=2.0.4,!=2.1.2,!=2.1.6,>=2.0.1 in /opt/conda/lib/python3.6/site-packages (from matplotlib>=1.4->pandas_profiling)
Installing collected packages: pandas-profiling
Successfully installed pandas-profiling-1.4.1
You are using pip version 9.0.1, however version 10.0.1 is available.
You should consider upgrading via the 'pip install --upgrade pip' command.
Collecting xgboost
  Downloading https://files.pythonhosted.org/packages/4b/c4/57e246bc99e45c048f9805f2773e7369f0d30896d19fa089fa1794c7b246/xgboost-0.71.tar.gz (494kB)
    100% |████████████████████████████████| 501kB 1.2MB/s ta 0:00:01
Requirement already satisfied: numpy in /opt/conda/lib/python3.6/site-packages (from xgboost)
Requirement already satisfied: scipy in /opt/conda/lib/python3.6/site-packages (from xgboost)
Building wheels for collected packages: xgboost
  Running setup.py bdist_wheel for xgboost ... done
  Stored in directory: /root/.cache/pip/wheels/4e/6d/1d/0bc23240225fe411315d8abb5d4521b9ff002493ff77515ccc
Successfully built xgboost
Installing collected packages: xgboost
Successfully installed xgboost-0.71
You are using pip version 9.0.1, however version 10.0.1 is available.
You should consider upgrading via the 'pip install --upgrade pip' command.

#import libraries
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import pandas_profiling
import seaborn as sns


from sklearn.linear_model import LinearRegression
from sklearn.ensemble import RandomForestRegressor
from sklearn.svm import SVR
from xgboost import XGBRegressor


from sklearn.preprocessing import Imputer,OneHotEncoder,LabelEncoder,StandardScaler
from sklearn.cross_validation import train_test_split
from sklearn.model_selection import cross_val_score
from sklearn.decomposition import PCA

Approach #1

This is going to be a very simple and naive approach of what we already know.
It is also useful to compare it with other later approaches.
We will proceed in the following steps:

Descriptive Statistics of the data
Data Cleaning and Feature Selection
Testing different models and comparing

Descriptive Statistics of the data

We begin by visualising the data to get some more feel of it.
We will import the data first.

def import_data():
    train=pd.read_csv("Challenge Data/train.csv")
    test=pd.read_csv("Challenge Data/test.csv")
    return train,test

train,test=import_data()

Let's get a sense of the data by looking at some rows.

train.head(5)

test.head(5)

We will use a very cool trick to visualise the input data taught from here.
Please note that clicking on "Toggle details" doesn't work unless you use an ipynb file.

pandas_profiling.ProfileReport(train)

We can see that we have some missing data,both for numerical and categorical data.
For categorical data,missing data refers to the fact that the house doesn't have the feature. For example, we have a lot of missing values in the feature PoolQC. This means that there are lots of houses that have a pool quality of "NaN",meaning they don't have a pool.
What we want to solve is the missing values for numerical features such as Lot Frontage.
But before we do so,let us first separate the features from the output and join both train and test dataframes so that we can preprocess their data together.

Data Cleaning and Feature Selection

#get the IDs
train_ID=train.Id
test_ID=test.Id

#get sales price
y=train.SalePrice

#drop the non-feature columns
train=train.drop(['Id','SalePrice'],axis=1)
test=test.drop('Id',axis=1)

We now append the test to the train set,so that we can preprocess their data together,but first we save a variable so that we can split the data then.

index=train.shape[0]

#combine both dataframes
combined=train.append(test)

Let us now deal with the missing values for numerical features.

#this is the number of elements for each feature
combined.describe().iloc[0,:]

MSSubClass       1460.0
LotFrontage      1201.0
LotArea          1460.0
OverallQual      1460.0
OverallCond      1460.0
YearBuilt        1460.0
YearRemodAdd     1460.0
MasVnrArea       1452.0
BsmtFinSF1       1460.0
BsmtFinSF2       1460.0
BsmtUnfSF        1460.0
TotalBsmtSF      1460.0
1stFlrSF         1460.0
2ndFlrSF         1460.0
LowQualFinSF     1460.0
GrLivArea        1460.0
BsmtFullBath     1460.0
BsmtHalfBath     1460.0
FullBath         1460.0
HalfBath         1460.0
BedroomAbvGr     1460.0
KitchenAbvGr     1460.0
TotRmsAbvGrd     1460.0
Fireplaces       1460.0
GarageYrBlt      1379.0
GarageCars       1460.0
GarageArea       1460.0
WoodDeckSF       1460.0
OpenPorchSF      1460.0
EnclosedPorch    1460.0
3SsnPorch        1460.0
ScreenPorch      1460.0
PoolArea         1460.0
MiscVal          1460.0
MoSold           1460.0
YrSold           1460.0
Name: count, dtype: float64

We can see that we have missing values in the features: LotFrontage,MasVnrArea,and GarageYrBlt.
Let us check each feature separatly and figure out how to impute the missing values.

combined.LotFrontage.describe()

count    1201.000000
mean       70.049958
std        24.284752
min        21.000000
25%        59.000000
50%        69.000000
75%        80.000000
max       313.000000
Name: LotFrontage, dtype: float64

We have a few undefined values in the columns. We are going to fill them by the mean of the column.
(filling by the mean is not the best approach,but since we are implementing a naive approach,we will fill by the mean as it is the one of the most common ways for imputing missing data)

combined.LotFrontage.fillna(combined.LotFrontage.mean(),inplace=True)

combined.LotFrontage.describe()

count    1460.000000
mean       70.049958
std        22.024023
min        21.000000
25%        60.000000
50%        70.049958
75%        79.000000
max       313.000000
Name: LotFrontage, dtype: float64

Let's check for MasVnArea.

combined.MasVnrArea.describe()

count    1452.000000
mean      103.685262
std       181.066207
min         0.000000
25%         0.000000
50%         0.000000
75%       166.000000
max      1600.000000
Name: MasVnrArea, dtype: float64

We also have a few missing values that we are going to fill with the mean of the column.

combined.MasVnrArea.fillna(combined.MasVnrArea.mean(),inplace=True)

Let's look at the last feature.

combined.GarageYrBlt.describe()

count    1379.000000
mean     1978.506164
std        24.689725
min      1900.000000
25%      1961.000000
50%      1980.000000
75%      2002.000000
max      2010.000000
Name: GarageYrBlt, dtype: float64

We also have a few missing values that we are going to fill with the mean of the column.

combined.GarageYrBlt.fillna(combined.GarageYrBlt.mean(),inplace=True)

combined.describe().iloc[0,:]

MSSubClass       1460.0
LotFrontage      1460.0
LotArea          1460.0
OverallQual      1460.0
OverallCond      1460.0
YearBuilt        1460.0
YearRemodAdd     1460.0
MasVnrArea       1460.0
BsmtFinSF1       1460.0
BsmtFinSF2       1460.0
BsmtUnfSF        1460.0
TotalBsmtSF      1460.0
1stFlrSF         1460.0
2ndFlrSF         1460.0
LowQualFinSF     1460.0
GrLivArea        1460.0
BsmtFullBath     1460.0
BsmtHalfBath     1460.0
FullBath         1460.0
HalfBath         1460.0
BedroomAbvGr     1460.0
KitchenAbvGr     1460.0
TotRmsAbvGrd     1460.0
Fireplaces       1460.0
GarageYrBlt      1460.0
GarageCars       1460.0
GarageArea       1460.0
WoodDeckSF       1460.0
OpenPorchSF      1460.0
EnclosedPorch    1460.0
3SsnPorch        1460.0
ScreenPorch      1460.0
PoolArea         1460.0
MiscVal          1460.0
MoSold           1460.0
YrSold           1460.0
Name: count, dtype: float64

We can see that we managed fix the issue of the missing values.
Let's now check the correlation between the features and the Sale Price feature.

combined['SalePrice']=y

combined.corr().SalePrice.sort_values(ascending=False)

SalePrice        1.000000
OverallQual      0.654385
GrLivArea        0.590617
GarageCars       0.524004
GarageArea       0.513377
TotalBsmtSF      0.501647
1stFlrSF         0.490801
FullBath         0.461570
TotRmsAbvGrd     0.448564
YearBuilt        0.432041
YearRemodAdd     0.427852
MasVnrArea       0.403802
GarageYrBlt      0.386839
Fireplaces       0.372401
BsmtFinSF1       0.333616
LotFrontage      0.290158
2ndFlrSF         0.278702
WoodDeckSF       0.270388
OpenPorchSF      0.258385
HalfBath         0.234466
LotArea          0.233986
BsmtFullBath     0.198040
BsmtUnfSF        0.146941
BedroomAbvGr     0.141105
ScreenPorch      0.104810
PoolArea         0.046896
MoSold           0.041965
BsmtFinSF2       0.018506
3SsnPorch        0.012448
BsmtHalfBath    -0.014629
LowQualFinSF    -0.021420
YrSold          -0.036494
MiscVal         -0.050300
OverallCond     -0.056972
MSSubClass      -0.059659
EnclosedPorch   -0.116645
KitchenAbvGr    -0.129535
Name: SalePrice, dtype: float64

The features PoolArea, MoSold, BsmtFinSF2, 3SsnPorch, BsmtHalfBath, LowQualFinSF, YrSold,MiscVal, OverallCond, and MSSubClass have very low correlations,so we are going to drop them.

combined.drop(['PoolArea','MoSold','BsmtFinSF2','3SsnPorch','BsmtHalfBath','LowQualFinSF','YrSold','MiscVal','OverallCond','MSSubClass','SalePrice'],axis=1,inplace=True)

Now time for the categorical features. We are going to one-hot encode them,but before that we need to fill the missing values by "None".

combined.Alley.isna().sum()

1369

combined.isna().sum().sort_values(ascending=False)

PoolQC           1453
MiscFeature      1406
Alley            1369
Fence            1179
FireplaceQu       690
GarageFinish       81
GarageType         81
GarageCond         81
GarageQual         81
BsmtFinType2       38
BsmtExposure       38
BsmtQual           37
BsmtCond           37
BsmtFinType1       37
MasVnrType          8
Electrical          1
Exterior1st         0
Exterior2nd         0
RoofMatl            0
ExterCond           0
MasVnrArea          0
ExterQual           0
YearRemodAdd        0
Foundation          0
RoofStyle           0
SaleCondition       0
YearBuilt           0
OverallQual         0
HouseStyle          0
BldgType            0
                 ... 
LotFrontage         0
Condition1          0
TotalBsmtSF         0
BsmtFinSF1          0
KitchenQual         0
ScreenPorch         0
EnclosedPorch       0
OpenPorchSF         0
WoodDeckSF          0
PavedDrive          0
GarageArea          0
GarageCars          0
GarageYrBlt         0
Fireplaces          0
Functional          0
TotRmsAbvGrd        0
KitchenAbvGr        0
BsmtUnfSF           0
BedroomAbvGr        0
HalfBath            0
FullBath            0
BsmtFullBath        0
GrLivArea           0
2ndFlrSF            0
1stFlrSF            0
CentralAir          0
HeatingQC           0
Heating             0
SaleType            0
MSZoning            0
Length: 69, dtype: int64

combined.fillna("None",inplace=True)

#hot encode the data
combined=pd.get_dummies(combined)

Let's check if we have any missing values.

combined.isnull().sum().sort_values(ascending=False)

SaleCondition_Partial    0
HouseStyle_SLvl          0
Condition2_PosN          0
Condition2_RRAe          0
Condition2_RRAn          0
Condition2_RRNn          0
BldgType_1Fam            0
BldgType_2fmCon          0
BldgType_Duplex          0
BldgType_Twnhs           0
BldgType_TwnhsE          0
HouseStyle_1.5Fin        0
HouseStyle_1.5Unf        0
HouseStyle_1Story        0
HouseStyle_2.5Fin        0
HouseStyle_2.5Unf        0
HouseStyle_2Story        0
Condition2_PosA          0
Condition2_Norm          0
Condition2_Feedr         0
Condition1_Feedr         0
Neighborhood_SawyerW     0
Neighborhood_Somerst     0
Neighborhood_StoneBr     0
Neighborhood_Timber      0
Neighborhood_Veenker     0
Condition1_Artery        0
Condition1_Norm          0
Condition2_Artery        0
Condition1_PosA          0
                        ..
Heating_OthW             0
Heating_Wall             0
HeatingQC_Ex             0
HeatingQC_Gd             0
BsmtFinType2_Unf         0
HeatingQC_Po             0
HeatingQC_TA             0
CentralAir_N             0
CentralAir_Y             0
Electrical_FuseA         0
Electrical_FuseF         0
Heating_Floor            0
BsmtFinType2_Rec         0
BsmtExposure_Gd          0
BsmtFinType1_LwQ         0
BsmtExposure_Mn          0
BsmtExposure_No          0
BsmtExposure_None        0
BsmtFinType1_ALQ         0
BsmtFinType1_BLQ         0
BsmtFinType1_GLQ         0
BsmtFinType1_None        0
BsmtFinType2_None        0
BsmtFinType1_Rec         0
BsmtFinType1_Unf         0
BsmtFinType2_ALQ         0
BsmtFinType2_BLQ         0
BsmtFinType2_GLQ         0
BsmtFinType2_LwQ         0
LotFrontage              0
Length: 293, dtype: int64

Let's now get the modified train and test sets.

mtrain=combined.iloc[:index,:]

mtest=combined.iloc[index:,:]

Testing Different Models and Comparing

We will train different models but we must first assess them.
We will apply cross fold validation in order to do so.

#this function takes an input as a model and prints the average rmse of k-fold cross validation
def model_score(model,mtrain,y,cv=50):
    #mse are by default negative in sklearn
    mse=-cross_val_score(model,mtrain,y,cv=50,scoring='neg_mean_squared_error')
    #take square root
    rmse=np.sqrt(mse)
    #average
    average_rmse=np.mean(rmse)
    print("The average RMSE is: %f." %rmse[0])

Now time to test some models.
We first initiate the model, then we fit the data to it, and finally we use model_score to assess the model.

#Linear Regressor
regressor=LinearRegression()
regressor.fit(mtrain,y)
model_score(regressor,mtrain,y)

The average RMSE is: 28940.975976.

#random forest
rf=RandomForestRegressor(n_estimators=100)
rf.fit(mtrain,y)
model_score(rf,mtrain,y)

The average RMSE is: 24086.086256.

#XGBoost
xgb=XGBRegressor()
xgb.fit(mtrain,y)
model_score(xgb,mtrain,y)

The average RMSE is: 22355.793436.

#SVR
svr=SVR(kernel='rbf')
svr.fit(mtrain,y)
model_score(svr,mtrain,y)

The average RMSE is: 72607.925589.

#we will the same models but we will apply PCA to the data
pca=PCA(n_components=100)
pca_mtrain=pca.fit_transform(mtrain)

regressor.fit(pca_mtrain,y)
model_score(regressor,pca_mtrain,y)

The average RMSE is: 28435.432337.

rf.fit(pca_mtrain,y)
model_score(rf,pca_mtrain,y)

The average RMSE is: 23866.225774.

xgb.fit(pca_mtrain,y)
model_score(xgb,pca_mtrain,y)

The average RMSE is: 25108.068674.

svr=SVR(kernel='rbf')
svr.fit(pca_mtrain,y)
model_score(svr,mtrain,y)

The average RMSE is: 72607.925589.

Remarks and Notes

Clearly this approach isn't the best as many things we did were not justified scientifically.
Even if you use the most sophisticated model, if you input garbage,then you will also output garbage.
A lot of work needs to be done on data cleaning and feature selection

All this will be tackled in the second approach.

Approach #2

The approach done before was a bit naive and lacks any statistical reasons behind it.
We fix this here in this approach and hope to get a lower rmse.

We first begin by reimporting the data again.

def import_data():
    train=pd.read_csv("Challenge Data/train.csv")
    test=pd.read_csv("Challenge Data/test.csv")
    return train,test


train,test=import_data()

Let's separate the IDs and Sale Price

#get IDs
train_ID=train.Id
test_ID=test.Id

#get SalePrices
y=train.SalePrice

train.shape

(1200, 81)

#drop the unwanted columns
train.drop(['Id','SalePrice'],axis=1,inplace=True)
test.drop('Id',axis=1,inplace=True)

train.shape

(1200, 79)

Let's see how much missing values we have.

#total number of rows of dataset
n_rows=train.shape[0]+test.shape[0]

#proportions of number of rows that have missing data to total number of rows 
test_miss=test.isna().sum()/n_rows
train_miss=train.isna().sum()/n_rows

plt.figure(figsize=(13,7))
plt.title("Percentage of missing values in training dataset")
train_miss[train_miss>0].sort_values(ascending=True).plot(kind="barh")
plt.xlabel("Percentage")
plt.show()

plt.figure(figsize=(13,7))
plt.title("Percentage of missing values in test dataset")
test_miss[test_miss>0].sort_values(ascending=True).plot(kind="barh")
plt.xlabel("Percentage")
plt.show()

In addition to the data visualisation done in Approach#1, we will explore the data more.
Let's start by plotting the Sale Price

plt.figure(figsize=(13,7))
sns.distplot(y)
plt.title("Distribution Plot of SalePrice")

Text(0.5,1,'Distribution Plot of SalePrice')

We see here that the SalePrice is skewed and has a peak value. This mean more people tend to buy houses which are not expensive and a few people buy expensive houses.
Let's check its Kurtosis and Skewness.
To simply put it, a higher kurtosis means that the distribution has fat tails compared to a standard normal distribution.
A Positive skeweness means that the distribution is positively skewed,i.e., the tail is located to the right.

print("Kurtosis: %f" %y.kurt())
print("Skewness: %f" %y.skew())

Kurtosis: 7.033907
Skewness: 1.967215

Apparently, many people tend to take the log of right skewed depenedent variables and work with it and then use the exponentiation function to get the result.
We will do so and compare both ways to see which is better.

plt.figure(figsize=(13,7))
log_y=np.log(y)
sns.distplot(log_y)
plt.show()
print("Kurtosis: %f" %log_y.kurt())
print("Skewness: %f" %log_y.skew())

Kurtosis: 0.888850
Skewness: 0.132714

We can see that the data has become more uniform.
We will now study numerical and categorical features.
Let's write a function that extracts the numerical and categorical features.

def get_features(train):
    #select all numerical feature names
    numerical_features=train.select_dtypes(include=np.number).columns

    #select all categorical feature names
    categorical_features=train.select_dtypes(exclude=np.number).columns
    
    return numerical_features,categorical_features

numerical_features,categorical_features=get_features(train)

We begin with numerical features. We will study each one by plotting distribution of each and a scatter plot against the SalePrice.
Let's begin with the distribution plot.

train.hist(figsize=(25,20))
plt.show()

A few remarks:

We have some features such as Garage Cars that seem to have a certain dominant value over the other values. A suggestion here might be to replace this dominant values by a "1" and all others by a "0".
It might be a good idea to sum all features related to a house area into one feature that denotes the total house area.
The YrSold variable seems to have nearly the same number of houses sold in all years except for 2010. That was probably due to the financial crisis that hit the world. We might encode it as a cat egorical variable with a value of "1" for years where sold houses are in 2010 and 0 for years which houses were not sold in 2010.
The GarageYrBuilt and the YearBuilt features are somehow similar. Usually,one builds the house along with the garage. It seems reasonable here to drop one of these variables as they are redundant.
The MiscVal feature seems useless as most of the data is 0 and it would be a good idea to drop the feature.
A lot of houses don't seem to have pools
It might be a good idea to sum features corresponding to the number of rooms,number of kitchens,number of bedrooms,and number of baths into one feature

We will tackle these in the Feature Engineering and Selection part.

Now let's plot scatter plots between each numerical features and the SalePrice variable.
As garage cars increase,garage area increases.

#dictionary contraing features name as keys and correlation with sale price as value
corr_dict=train[numerical_features].corrwith(y).sort_values(ascending=False).to_dict()

#graphs are sorted by alphabetical order of the feature names
i=1
plt.figure(figsize=(55, 60))
for k in sorted(corr_dict):
    plt.subplot(6,6,i)
    plt.scatter(train[k],y)
    plt.title("Correlation: %f" %corr_dict[k],fontsize=30)
    plt.xlabel(k,fontsize=25)
    plt.ylabel("Sale Prices",fontsize=25)
    i+=1

A few remarks:

There are some outliers in features such as GrLivArea and TotalBsmtSF. We should consider removing them.
There are features such as BsmtFinSF2 that have near zero correlation. In the previous approach, we removed such low corrleations. In this approach, we will keep a few of them and try to generate some useful features that might help us.

To finish visualisation of numerical features,let us plot a heat map to see correlations between the features.

#tryong colors for heat map
sns.palplot(sns.color_palette("RdBu_r", 6))

plt.subplots(figsize=(15,10))
ax=sns.heatmap(train.corr(),cmap=sns.color_palette("RdBu_r", 5),linewidths=.5,vmin=-1,vmax=1)
cbar = ax.collections[0].colorbar
cbar.set_ticks(np.arange(-1,1.1,0.2))
plt.title("Heatmap of Numerical Features")
plt.show()

We see that we have a lot of correlations between -0.2 and 0.2.
There also exists some nearly perfect corrleations between some features.
We will now study categorical features and then try to find best features to match our data.

We now study categorical features.

i=1
plt.figure(figsize=(80, 350))
for cat in sorted(categorical_features):
    plt.subplot(22,2,i)
    
    #this is just for the xlabels to not overlap
    if(cat=="Neighborhood"):
        plt.xticks(rotation=30)
        
    
    #stripplot generalizes scatter to categorical variables
    sns.stripplot(train[cat],y)
    plt.xlabel(cat,fontsize=30)
    plt.ylabel("Sale Prices",fontsize=30)
    plt.xticks(fontsize=27)
    plt.yticks(fontsize=27)
    i+=1

Obviously,we have missing values by looking at the feature PoolQC. We will deal with those in the next section.
Most of the data in the feature Utilities is dominated by a single value and we could consider removing it. Usually,people buy a house equiped with all the utilities.(Electricity, Water, Gas, and Sewer)

We will now impute the missing data.

Data Cleaning

Let's see how much missing values we have for numerical features in the train and test dataset.

miss_train=train.isna().sum().sort_values(ascending=False)
miss_train[miss_train>0]

PoolQC          1196
MiscFeature     1153
Alley           1125
Fence            973
FireplaceQu      564
LotFrontage      210
GarageType        67
GarageFinish      67
GarageQual        67
GarageCond        67
GarageYrBlt       67
BsmtExposure      33
BsmtFinType2      33
BsmtCond          32
BsmtFinType1      32
BsmtQual          32
MasVnrArea         6
MasVnrType         6
dtype: int64

miss_test=test.isna().sum().sort_values(ascending=False)
miss_test[miss_test>0]

PoolQC          257
MiscFeature     253
Alley           244
Fence           206
FireplaceQu     126
LotFrontage      49
GarageCond       14
GarageType       14
GarageYrBlt      14
GarageFinish     14
GarageQual       14
BsmtFinType1      5
BsmtExposure      5
BsmtCond          5
BsmtQual          5
BsmtFinType2      5
MasVnrArea        2
MasVnrType        2
Electrical        1
dtype: int64

Let's combine the train and test set so that we can clean thier data together.

#used to slice test from train
index=train.shape[0]

index

1200

full=pd.concat([train,test],ignore_index=True)

Let's start with numerical features.

miss_full=full[numerical_features].isna().sum().sort_values(ascending=False)
miss_full=miss_full[miss_full>0]
miss_full

LotFrontage    259
GarageYrBlt     81
MasVnrArea       8
dtype: int64

miss_full_cols=miss_full.index

pandas_profiling.ProfileReport(full[miss_full_cols])

For the feature GarageYrBlt, we have some missing values. What we will do now is fill the missing values by zero since this means that these houses have no garage. We will experiment later on with other approaches such as removing the rows or filling them with the mode.
But before that,let's look more at the GarageYrBlt feature. We have that the minimum value for the feature "GarageYrBuilt" is 1900 whereas the minimum value for "YearBuilt" is 1872. That seems weird. Investigating more in the feature we observe that:

We have some rows where the garage was built many years after the house was built. This seems fine as the house might have been passed along many generations and one of those generations decided to buy a car and thus build a garage.
We have some rows where the garage was built before the house was built. Although this might seem strange,but if we started building a house and a garage,the garage clearly will be finished earlier than the house, so the year where the garage was finished must be less than that when the house was finished.

full[full.YearBuilt==1880][['YearBuilt','GarageYrBlt']]

full[full.GarageYrBlt<full.YearBuilt][['YearBuilt','GarageYrBlt']]

plt.figure(figsize=(13,9))
full.GarageYrBlt.hist()
plt.title("Distribution Plot of GarageYrBlt")
plt.show()

#fill missing values by 0
full.GarageYrBlt.fillna(0,inplace=True)

For the feature LotFrontage,we are going to fill the missing values with the mode.

full.LotFrontage.fillna(full.LotFrontage.mode()[0],inplace=True)

For the feature MasVnrArea, we will fill the missing values by zero as it is the most common element.

full.MasVnrArea.fillna(0,inplace=True)

Before tackling categorical features,let us check that we don't have any missing values.

full[numerical_features].isnull().sum().sort_values(ascending=False)>0

YrSold           False
MoSold           False
GrLivArea        False
LowQualFinSF     False
2ndFlrSF         False
1stFlrSF         False
TotalBsmtSF      False
BsmtUnfSF        False
BsmtFinSF2       False
BsmtFinSF1       False
MasVnrArea       False
YearRemodAdd     False
YearBuilt        False
OverallCond      False
OverallQual      False
LotArea          False
LotFrontage      False
BsmtFullBath     False
BsmtHalfBath     False
FullBath         False
WoodDeckSF       False
MiscVal          False
PoolArea         False
ScreenPorch      False
3SsnPorch        False
EnclosedPorch    False
OpenPorchSF      False
GarageArea       False
HalfBath         False
GarageCars       False
GarageYrBlt      False
Fireplaces       False
TotRmsAbvGrd     False
KitchenAbvGr     False
BedroomAbvGr     False
MSSubClass       False
dtype: bool

Now it's time for the categorical variables. What we will do now is to replace every "nan" with "None" because when applying dummy enconding "None" is considered as a column while "nan" is not.

#number of missing values for categorical features
miss_full_cat=full[categorical_features].isna().sum().sort_values(ascending=False)
miss_full_cat=miss_full_cat[miss_full_cat>0]
miss_full_cat

PoolQC          1453
MiscFeature     1406
Alley           1369
Fence           1179
FireplaceQu      690
GarageCond        81
GarageQual        81
GarageFinish      81
GarageType        81
BsmtFinType2      38
BsmtExposure      38
BsmtFinType1      37
BsmtQual          37
BsmtCond          37
MasVnrType         8
Electrical         1
dtype: int64

Consulting the "Data description" file given,most categorical features that contain missing values correspond to not having that feature. For example,missing values for the feature "PoolQC" means that the house doesn't contain a pool.
However,there is only one categorical feature that contains "nan" values and it shouldn't which is the "Electrical" feature
Let's deal with it.

full["Electrical"].describe()

count      1459
unique        5
top       SBrkr
freq       1334
Name: Electrical, dtype: object

miss_full_cat['Electrical']

1

print("Unique Values for the feature Electrical:", full["Electrical"].unique())
print("Number of missing values for the feature Electrical:", miss_full_cat["Electrical"])
print("\n")

Unique Values for the feature Electrical: ['SBrkr' 'FuseF' 'FuseA' 'FuseP' 'Mix' nan]
Number of missing values for the feature Electrical: 1

plt.figure(figsize=(13,9))
full.Electrical.hist()
plt.title("Distribution plot of Electrical")
plt.show()

full[full.Electrical.isnull()>0]

Since the number of missing data is small,it doesn't hurt fill it with the most common element "SBrkr".

full.Electrical.fillna('SBrkr',inplace=True)

Now having finished all this we can proceed for feature engineering.

Feature Engineering and Selection

Let's start with the numerical features.
Let's replace the features that have to do with house area with one feature corresponding to the totalSF.

Let's have a look at the the features related to basement area.
We have the features BsmtFinSF1,BsmtFinSF2, BsmtUnfSF, and TotalBsmtSF.
It wouldn't be a surprise if the sum of the first three features is equal to data in TotalBsmtSF.
Let's see!

#sum all basement area
bsmtsum=full.BsmtFinSF1+full.BsmtFinSF2+full.BsmtUnfSF

#compare to total sum
bsmtsum.equals(full.TotalBsmtSF)

True

This means that we can drop the 3 features and keep only the TotalBsmtSF feature!

full.shape

(1460, 79)

full.drop('BsmtFinSF1',axis=1,inplace=True)
full.drop('BsmtFinSF2',axis=1,inplace=True)
full.drop('BsmtUnfSF',axis=1,inplace=True)

Now let's look at the area of the house. The feature GrLivArea is the area above ground.
Let's try to sum some features and see.

(full['1stFlrSF']+full['2ndFlrSF']+full.LowQualFinSF).equals(full.GrLivArea)

True

It seems that the sum of 1stFlrSF, 2ndFlrSF, and LowQualFinSF equals the GrLivArea feature.
We can keep only the GrLivArea and drop the others.

full.drop('LowQualFinSF',axis=1,inplace=True)
full.drop('1stFlrSF',axis=1,inplace=True)
full.drop('2ndFlrSF',axis=1,inplace=True)

Let's look at the porches.

(full['3SsnPorch']+full.EnclosedPorch+full.OpenPorchSF+full.ScreenPorch).iloc[:index].corr(y)

0.19782953308159243

We still have these raw features left: 3SsnPorch,EnclosedPorch,GarageArea,LotArea,MasVnrArea,OpenPorchSF,PoolArea,ScreenPorch,WoodDeckSF. We will leave them for now.
Let's drop useless features such MiscVal.

full.drop('MiscVal',axis=1,inplace=True)

Now off to categorical features.
Let's plot histograms of each feature.

i=1
plt.figure(figsize=(50, 250))
for cat in sorted(categorical_features):
    plt.subplot(22,2,i)

    
    #this is just for the xlabels to not overlap
    if(cat in ["Neighborhood","Exterior1st","Exterior2nd"]):
        plt.xticks(rotation=30)    
    

    full[cat].hist()
    plt.xticks(fontsize=20)
    plt.yticks(fontsize=20)
    plt.title("{}".format(cat),size=20)

    #plt.show()
    i+=1

Exterior1st looks much like Exterior2nd so it is safe to remove one of them.

full.drop('Exterior1st',axis=1,inplace=True)

We can remove Utilities because nearly all values are the same and that doens't help us much.

full.drop('Utilities',axis=1,inplace=True)

full.isnull().sum().sort_values(ascending=False)

PoolQC           1453
MiscFeature      1406
Alley            1369
Fence            1179
FireplaceQu       690
GarageType         81
GarageCond         81
GarageFinish       81
GarageQual         81
BsmtFinType2       38
BsmtExposure       38
BsmtCond           37
BsmtQual           37
BsmtFinType1       37
MasVnrType          8
LandSlope           0
LotConfig           0
MSZoning            0
LotFrontage         0
LotArea             0
Foundation          0
ExterCond           0
ExterQual           0
MasVnrArea          0
Exterior2nd         0
RoofMatl            0
RoofStyle           0
YearRemodAdd        0
YearBuilt           0
OverallCond         0
                 ... 
HeatingQC           0
TotalBsmtSF         0
Fireplaces          0
YrSold              0
MoSold              0
PoolArea            0
ScreenPorch         0
3SsnPorch           0
EnclosedPorch       0
OpenPorchSF         0
WoodDeckSF          0
PavedDrive          0
GarageArea          0
GarageCars          0
GarageYrBlt         0
Functional          0
Heating             0
TotRmsAbvGrd        0
KitchenQual         0
KitchenAbvGr        0
BedroomAbvGr        0
HalfBath            0
FullBath            0
BsmtHalfBath        0
BsmtFullBath        0
GrLivArea           0
Electrical          0
CentralAir          0
SaleType            0
MSSubClass          0
Length: 70, dtype: int64

Now let's fill the missing values by "None".

full.fillna("None",inplace=True)

What we will do now is hot encode the data. We will try with removing one columns of each dummy variable to avoid the dummy trap just for experimenting

full_wd=full.copy() #wd=with droping first dummy variable column
full_wod=full.copy()#wod=without dropping first dummy variable column

#dummy encode
full_wd=pd.get_dummies(full_wd)
full_wod=pd.get_dummies(full_wod,drop_first=True)

#split to train and test
train_wd=full_wd.iloc[:index,:]
test_wd=full_wd.iloc[index:,:]

#split to train and test
train_wod=full_wod.iloc[:index,:]
test_wod=full_wod.iloc[index:,:]

Training Time

from sklearn.cross_validation import cross_val_score
def model_score(model,mtrain,y,cv=50):
    #mse are by default negative in sklearn
    mse=-cross_val_score(model,mtrain,y,cv=cv,scoring='neg_mean_squared_error')
    rmse=np.sqrt(mse)
    average_rmse=np.mean(rmse)
    std=np.std(rmse)
    #print(rmse)
    print("The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

We will try several models.

#linear regressor: train_wd+log_y
lm1=LinearRegression()
lm1.fit(train_wd,log_y)
model_score(lm1,train_wd,log_y)

The average RMSE is: 0.137537 with a standard deviation of 0.057180.

#linear regressor: train_wod+log_y
lm2=LinearRegression()
lm2.fit(train_wod,log_y)
model_score(lm2,train_wod,log_y)

The average RMSE is: 0.133158 with a standard deviation of 0.051099.

Logging output variable and removing dummy variable seems to work best.

#random forest: train_wod+log_y
rlf=RandomForestRegressor()
rlf.fit(train_wod,log_y)
model_score(rlf,train_wod,log_y)

The average RMSE is: 0.149585 with a standard deviation of 0.036118.

#random forest: train_wd+log_y
rlf2=RandomForestRegressor()
rlf2.fit(train_wd,log_y)
model_score(rlf2,train_wd,log_y)

The average RMSE is: 0.149487 with a standard deviation of 0.038988.

We will now use XGBoost which is a kind of boosted trees.

#XGBoost: train_wod+log_y
xgb=XGBRegressor()
xgb.fit(train_wod,log_y)
model_score(xgb,train_wod,log_y)

The average RMSE is: 0.129749 with a standard deviation of 0.038742.

#XGBoost: train_wdd+log_y
xgb=XGBRegressor()
xgb.fit(train_wd,log_y)
model_score(xgb,train_wd,log_y)

The average RMSE is: 0.129172 with a standard deviation of 0.038268.

Our approach this time produced better results.
However we seek an even less RMSE.
Aside from the models,the decisive factor for a lower RMSE is choosing the right features and we didn't do a logical plan to choose them in this approach.We fix this in the following approach.

Approach #3

Let us first execute the essential starting code.

train,test=import_data()

#get IDs
train_ID=train.Id
test_ID=test.Id

#get SalePrices
y=train.SalePrice

#drop the unwanted columns
train.drop(['Id','SalePrice'],axis=1,inplace=True)
test.drop('Id',axis=1,inplace=True)

#log the target
log_y=np.log(y)

numerical_features,categorical_features=get_features(train)

#combine train and test data
full=pd.concat([train,test],ignore_index=True)

We did some data visualisation in previous approaches and thus won't do it here.
We will directly work with data cleaning and feature selection. Let's look again at the missing values for numerical and categorcal features.
Let's start with numerical features.

num_miss=full[numerical_features].isna().sum().sort_values(ascending=False)
num_miss[num_miss>0]

LotFrontage    259
GarageYrBlt     81
MasVnrArea       8
dtype: int64

We have only 3 numerical features that have missing values.
Let's look at each one separately.
LotFrontage is the linear feet of street connected to property.
Usually, houses in the same neighborhood have their Lot Frontage values very close as seen below. So we will compute the average LotFrontage for each neighborhood and fill the missing values by the average.
For example,the average Lot Frontage for the Neighborhood Blmngtn is 47.142857,so we will fill each row that has Blmngtn as a Neighborhood and has a missing value for Lot Frontage by this average. Ofcourse,we are going to round this value as all values in LotFrontage are integers.

df=full[['LotFrontage','Neighborhood']]

neighb=df.Neighborhood.unique()

neighb

array(['CollgCr', 'Veenker', 'Crawfor', 'NoRidge', 'Mitchel', 'Somerst',
       'NWAmes', 'OldTown', 'BrkSide', 'Sawyer', 'NridgHt', 'NAmes',
       'SawyerW', 'IDOTRR', 'MeadowV', 'Edwards', 'Timber', 'Gilbert',
       'StoneBr', 'ClearCr', 'NPkVill', 'Blmngtn', 'BrDale', 'SWISU',
       'Blueste'], dtype=object)

for n in neighb:
    df[df.Neighborhood==n].plot.bar(figsize=(10,5),xticks=None,color = "skyblue")
    plt.tick_params(axis='x',which='both',bottom='off',top='off',labelbottom='off') 
    plt.title("Lot Frontages for %s" %n)
    plt.axhline(y=df[df.Neighborhood==n].mean()[0], color='r', linestyle='--',label='Mean')
    plt.ylabel("Lot Frontage")
    plt.legend()

As we can see, usually most of the houses have their Lot Frontage close to the mean,so it seems a good value to fill the missing value with.
Let's now fill the missing values with the average of each neighborhood.

full.LotFrontage=full.groupby("Neighborhood")["LotFrontage"].transform(lambda x: x.fillna(round(x.mean())))

#check for missing values
full.LotFrontage.isna().sum()

0

We now have to deal with the other 2 numerical variables.
Let's start with GarageYrBlt. These houses don't have garages so it seems logical that for these houses all features that have to do with garages are null.
Let's check.

garage_features=[n  for n in full.columns if 'Garage' in n]
garage_features

['GarageType',
 'GarageYrBlt',
 'GarageFinish',
 'GarageCars',
 'GarageArea',
 'GarageQual',
 'GarageCond']

full[garage_features].isna().sum()

GarageType      81
GarageYrBlt     81
GarageFinish    81
GarageCars       0
GarageArea       0
GarageQual      81
GarageCond      81
dtype: int64

All garage features have the same number of missing values except for Garage Cars and Garage Area.
Their values should sum to zero since if they don't have any garage than their garage cars and area must be zero. Also,their count should be 81.
Let us check.

full[full.GarageType.isna()].GarageCars.sum()

0

full[full.GarageCars==0].shape[0]

81

full[full.GarageType.isna()].GarageArea.sum()

0

full[full.GarageArea==0].shape[0]

81

Also, we need to check that the 81 missing values of the garage features correspond to the same entries.

full.GarageType.isna().equals(full.GarageYrBlt.isna())

True

full.GarageType.isna().equals(full.GarageCond.isna())

True

full.GarageType.isna().equals(full.GarageQual.isna())

True

Now let us discuss a procedure for filling the GarageYrBlt feature.
These houses don't have a garage so it would be useless to specify a year where the garage would be built.
If we fill these missing values with zero, we observe a big decrease in the correlation with the target feature.

full.GarageYrBlt.corr(y)

0.48979373237342116

full.GarageYrBlt.fillna(0).corr(y)

0.2551230905554551

However,looking at the graphs in a previous approach and consulting the correlation matrix, we find a large correlation with the feature YearBuilt.

full.GarageYrBlt.corr(full.YearBuilt)

0.825667484174342

Thus we are going to drop this feature.

full.drop('GarageYrBlt',axis=1,inplace=True)

We still have one numerical feature, MasVnrArea which is the masonry veneer area.
A missing value in this feature should mean that houses with these missing values don't have masonry veneers. Thus all features related to masonry veneer should be missing or have a value of zero.
Let's check.

mason_feat=[ind for ind in full.columns if 'Mas' in ind]

full[mason_feat].isna().sum()

MasVnrType    8
MasVnrArea    8
dtype: int64

It would be logical to fill houses with no masonry veneer area with 0.

full.MasVnrArea.fillna(0,inplace=True)

Let's check if we still have any missing values for numerical features aside from the GarageYrBlt feature.

numerical_features,categorical_features=get_features(full)

full[numerical_features].isna().sum().sort_values(ascending=False)

YrSold           0
BsmtFinSF1       0
LowQualFinSF     0
2ndFlrSF         0
1stFlrSF         0
TotalBsmtSF      0
BsmtUnfSF        0
BsmtFinSF2       0
MasVnrArea       0
BsmtFullBath     0
YearRemodAdd     0
YearBuilt        0
OverallCond      0
OverallQual      0
LotArea          0
LotFrontage      0
GrLivArea        0
BsmtHalfBath     0
MoSold           0
WoodDeckSF       0
MiscVal          0
PoolArea         0
ScreenPorch      0
3SsnPorch        0
EnclosedPorch    0
OpenPorchSF      0
GarageArea       0
FullBath         0
GarageCars       0
Fireplaces       0
TotRmsAbvGrd     0
KitchenAbvGr     0
BedroomAbvGr     0
HalfBath         0
MSSubClass       0
dtype: int64

Now it's time for categorical features.

full[categorical_features].isna().sum().sort_values(ascending=False)

PoolQC           1453
MiscFeature      1406
Alley            1369
Fence            1179
FireplaceQu       690
GarageCond         81
GarageQual         81
GarageFinish       81
GarageType         81
BsmtFinType2       38
BsmtExposure       38
BsmtFinType1       37
BsmtQual           37
BsmtCond           37
MasVnrType          8
Electrical          1
Condition2          0
Condition1          0
Neighborhood        0
LandSlope           0
BldgType            0
LandContour         0
LotConfig           0
Utilities           0
RoofStyle           0
LotShape            0
Street              0
HouseStyle          0
SaleCondition       0
RoofMatl            0
Exterior1st         0
Exterior2nd         0
ExterQual           0
ExterCond           0
Foundation          0
SaleType            0
Heating             0
HeatingQC           0
CentralAir          0
KitchenQual         0
Functional          0
PavedDrive          0
MSZoning            0
dtype: int64

Let's go through each feature beginning with PoolQC
Missing values in this feature means that the house has no pool,thus the pool area should be zero here.
Let's check.

full.PoolArea[full.PoolQC.isna()].sum()

0

Also we shouldn't have a Pool Area of zero when the PoolQC is not null.

full[full.PoolArea==0].equals(full[full.PoolQC.isna()])

True

We will now fill the missing values here with "No Pool".

full.PoolQC.fillna('No Pool',inplace=True)

Now let's see the feature MiscFeature.

plt.figure(figsize=(10,7))
full.MiscFeature.hist()
plt.xlabel("MiscFeature")
plt.ylabel("Count")
plt.show()

I will replace the missing values with "None".

full.MiscFeature.fillna("None",inplace=True)

plt.figure(figsize=(10,7))
full.MiscFeature.hist()
plt.xlabel("MiscFeature")
plt.ylabel("Count")
plt.show()

Now let's look at Alley.

plt.figure(figsize=(10,7))
full.Alley.hist()
plt.xlabel("Alley")
plt.ylabel("Count")
plt.show()

I will replace the missing value with "No Alley".

full.Alley.fillna("No Alley",inplace=True)

plt.figure(figsize=(10,7))
full.Alley.hist()
plt.xlabel("Alley")
plt.ylabel("Count")
plt.show()

The features Alley and MiscFeatures seem to have a lot of their data dominated by one value.
We might consider removing them later.

Now time for Fence.
We will fill the missing values by "No Fence".

plt.figure(figsize=(10,7))
full.Fence.fillna("No Fence",inplace=True)
full.Fence.hist()
plt.xlabel("Fence")
plt.ylabel("Count")
plt.show()

Now time for the feature FireplaceQu.

plt.figure(figsize=(10,7))
full.FireplaceQu.hist()
plt.xlabel("Fireplace Quality")
plt.ylabel("Count")
plt.show()

The houses with missing values for FireplaceQu should be identical to the rows of houses that have a zero in the feature Fireplaces.

full[full.FireplaceQu.isna()].equals(full[full.Fireplaces==0])

True

We will fill the missing values with "None".

full.FireplaceQu.fillna("NA",inplace=True)
plt.figure(figsize=(10,7))
full.FireplaceQu.hist()
plt.xlabel("Fireplace Quality")
plt.ylabel("Count")
plt.show()

We can replace the values here by numbers as "Ex" is greater than "Po" in terms of fire place quality.
We will do that in the feature engineering part.

Now it's time for the garage variables.
We will fill the missing values with "No Garage".

full.GarageCond.fillna('No Garage',inplace=True)
full.GarageQual.fillna('No Garage',inplace=True)
full.GarageFinish.fillna('No Garage',inplace=True)
full.GarageType.fillna('No Garage',inplace=True)

#just to check
full[full.GarageArea==0].equals(full[full.GarageCond=='No Garage'])

True

We can replace the values for the features GarageCond,GarageQual, and GarageFinish with numbers.
We will do so in the feature engineering part.

Now it's time for the basement variables.
We seem to have different numbers of missing values.
Let's check.

full.TotalBsmtSF.isna().sum()

0

bsmt_features=[ind for ind in full.columns if 'Bsmt' in ind]

full[bsmt_features].isna().sum().sort_values(ascending=False)

BsmtFinType2    38
BsmtExposure    38
BsmtFinType1    37
BsmtCond        37
BsmtQual        37
BsmtHalfBath     0
BsmtFullBath     0
TotalBsmtSF      0
BsmtUnfSF        0
BsmtFinSF2       0
BsmtFinSF1       0
dtype: int64

full[bsmt_features][full.BsmtFinType1.isna()].shape

(37, 11)

full[bsmt_features][full.BsmtFinType2.isna()].shape

(38, 11)

There is only one additional row.
Let's check it.

i=full[bsmt_features][full.BsmtFinType2.isna()].index.difference(full[bsmt_features][full.BsmtFinType1.isna()].index)[0]
i

332

pd.DataFrame(full[bsmt_features].iloc[i,:])

Here,we have a missing value for BsmtFinType2. Since,this house has a value in the feature BsmtFinSF2,it should have a value in BsmtFinType2.
Let's look at the feature BsmtFinType2.

full.BsmtFinType2.hist()
plt.show()

Since most of the values of this feature are "Unf",we think hat it might be a good value.
However,we can't do this here as we only use "Unf" when the value of BsmtFinSF2 is zero and this is not the case here as we have finished area.

full[['BsmtFinSF2','BsmtFinType2']][full.BsmtFinType2=='Unf']

[(ind,full.BsmtFinType2[full.BsmtFinType2==ind].count()) for ind in full.BsmtFinType2.unique()]

[('Unf', 1256),
 ('BLQ', 33),
 (nan, 0),
 ('ALQ', 19),
 ('Rec', 54),
 ('LwQ', 46),
 ('GLQ', 14)]

We see that the second most common value is Rec and so we will fill the missing value with this.

full.BsmtFinType2.iloc[i]='Rec'

Now we will see the additional missing value in BsmtExposure.

i=full[bsmt_features][full.BsmtExposure.isna()].index.difference(full[bsmt_features][full.BsmtFinType1.isna()].index)[0]
i

948

pd.DataFrame(full[bsmt_features].iloc[i,:])

Since this house has a basement,we can't fill the missing value with "No Basement".
Let's try to figure a way.

full.BsmtExposure.hist()
plt.show()

Since the most common value is "No" and there is no direct relation BsmtExposure and other basement features, we will fill the missing value by "No".

full.BsmtExposure.iloc[i]='No'

Now we should have the same missing values for the basement features.

full[bsmt_features].isna().sum()

BsmtQual        37
BsmtCond        37
BsmtExposure    37
BsmtFinType1    37
BsmtFinSF1       0
BsmtFinType2    37
BsmtFinSF2       0
BsmtUnfSF        0
TotalBsmtSF      0
BsmtFullBath     0
BsmtHalfBath     0
dtype: int64

#check if missing rows are the same for categorical features
full[full.BsmtQual.isna()].equals(full[full.BsmtFinType1.isna()])\
and full[full.BsmtQual.isna()].equals(full[full.BsmtFinType2.isna()])\
and full[full.BsmtQual.isna()].equals(full[full.BsmtCond.isna()]) \
and full[full.BsmtQual.isna()].equals(full[full.BsmtExposure.isna()])

True

#check for numerical variable
full[full.BsmtQual.isna()].equals(full[full.TotalBsmtSF==0])

True

Now we can fill the remaining values with "No Basement".

full.BsmtQual.fillna('No Basement',inplace=True)
full.BsmtCond.fillna('No Basement',inplace=True)
full.BsmtFinType1.fillna('No Basement',inplace=True)
full.BsmtExposure.fillna('No Basement',inplace=True)
full.BsmtFinType2.fillna('No Basement',inplace=True)

#check if we have any missing values for basement features
full[bsmt_features].isna().sum()

BsmtQual        0
BsmtCond        0
BsmtExposure    0
BsmtFinType1    0
BsmtFinSF1      0
BsmtFinType2    0
BsmtFinSF2      0
BsmtUnfSF       0
TotalBsmtSF     0
BsmtFullBath    0
BsmtHalfBath    0
dtype: int64

Let's now look at the last two features.
For the feature MasVnrType,we can directly fill the missing value with "No Masonry Veneer".

full.MasVnrType.fillna("No Masonry Veneer",inplace=True)

For the Electrical feature, we will fill the missing value with the most common value.

full.Electrical.hist()
plt.xlabel("Electrical")
plt.ylabel("Count")
plt.show()

full.Electrical.fillna('SBrkr',inplace=True)

#check to see if we have any missing values
full.isna().sum().sort_values(ascending=False)

SaleCondition    0
Foundation       0
RoofMatl         0
Exterior1st      0
Exterior2nd      0
MasVnrType       0
MasVnrArea       0
ExterQual        0
ExterCond        0
BsmtQual         0
YearRemodAdd     0
BsmtCond         0
BsmtExposure     0
BsmtFinType1     0
BsmtFinSF1       0
BsmtFinType2     0
BsmtFinSF2       0
BsmtUnfSF        0
RoofStyle        0
YearBuilt        0
SaleType         0
Utilities        0
MSZoning         0
LotFrontage      0
LotArea          0
Street           0
Alley            0
LotShape         0
LandContour      0
LotConfig        0
                ..
EnclosedPorch    0
ScreenPorch      0
CentralAir       0
PoolArea         0
PoolQC           0
Fence            0
MiscFeature      0
MiscVal          0
MoSold           0
YrSold           0
GarageCars       0
GarageFinish     0
GarageType       0
FireplaceQu      0
Electrical       0
1stFlrSF         0
2ndFlrSF         0
LowQualFinSF     0
GrLivArea        0
BsmtFullBath     0
BsmtHalfBath     0
FullBath         0
HalfBath         0
BedroomAbvGr     0
KitchenAbvGr     0
KitchenQual      0
TotRmsAbvGrd     0
Functional       0
Fireplaces       0
MSSubClass       0
Length: 78, dtype: int64

#we will keep a copy of the cleaned data for future use
full_cleaned=full.copy()

#full=full_cleaned.copy()

Feature Engineering

Let's begin with categorical features.
We will replace some categorcial features with numbers.
For example, let's look at the feature PoolQC.

full.PoolQC.unique()

array(['No Pool', 'Ex', 'Fa', 'Gd'], dtype=object)

We can see that 'Ex' is greater than all other values in terms of quality.
So we can replace the values by numbers with higher numbers meaning higher quality.

dict1={'No Pool':0,'Fa':1,'Gd':2,'Ex':3}

full.PoolQC=full.PoolQC.apply(lambda x: dict1[x])

Let's look at MiscFeatures.

full.MiscFeature.hist()
plt.xlabel("MiscFeatures")
plt.ylabel("Count")
plt.show()

There seems to be a dominant value here.
Let's try replacing the dominant value "None" by a zero and the others by a "1" and we will see its correlation with the target variable.
If it is too low,we will delete it.

a=full.MiscFeature.apply(lambda x: 0 if x=='None'  else 1)

a.corr(y)

-0.09204248012681689

delete=[] #list to store features to be deleted
delete.append("MiscFeature")

Let's look at Alley.

full.Alley.hist()
plt.xlabel("Alley")
plt.ylabel("Count")
plt.show()

There seems to be a dominant value here.
Let's try replacing the dominant value "None" by a zero and the others by a "1" and we will see its correlation with the target variable.
If it is too low,we will delete it.

a=full.Alley.apply(lambda x: 0 if x=='No Alley'  else 1)

a.corr(y)

-0.12220790878535244

delete.append('Alley')

Let's look at Fence.

full.Fence.hist()
plt.xlabel("Fence")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Let's look at the feature FireplaceQu.

full.FireplaceQu.hist()
plt.xlabel("Fireplace Quality")
plt.ylabel("Count")
plt.show()

We will replace these categories by numbers as we did before.

dict1={'NA':0,'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}

full.FireplaceQu=full.FireplaceQu.apply(lambda x: dict1[x])

Now let's look at the garage features and replace the necessary categories by numbers.

full.GarageCond.hist()
plt.xlabel("GarageCond")
plt.ylabel("Count")
plt.show()

#replace by numbers
dict1={'No Garage':0,'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}

full.GarageCond=full.GarageCond.apply(lambda x: dict1[x])

full.GarageQual.hist()
plt.xlabel("GarageQual")
plt.ylabel("Count")
plt.show()

full.GarageQual=full.GarageQual.apply(lambda x: dict1[x])

full.GarageFinish.hist()
plt.xlabel("GarageFinish")
plt.ylabel("Count")
plt.show()

#a finished garage cost more
dict1={'No Garage':0,'Unf':1,'RFn':2,'Fin':3}

full.GarageFinish=full.GarageFinish.apply(lambda x: dict1[x])

Now let's work with the basement features.

full.BsmtFinType2.unique()

array(['Unf', 'BLQ', 'No Basement', 'ALQ', 'Rec', 'LwQ', 'GLQ'],
      dtype=object)

#higher basement rating means higher price
dict1={'No Basement':0, 'Unf':1, 'LwQ':2, 'Rec':3, 'BLQ':4, 'ALQ':5, 'GLQ':6}

full.BsmtFinType1=full.BsmtFinType1.apply(lambda x: dict1[x])
full.BsmtFinType2=full.BsmtFinType2.apply(lambda x: dict1[x])

full.BsmtExposure.unique()

array(['No', 'Gd', 'Mn', 'Av', 'No Basement'], dtype=object)

dict1={'No Basement':0, 'No':1, 'Mn':2, 'Av':3, 'Gd':4}

full.BsmtExposure=full.BsmtExposure.apply(lambda x: dict1[x])

full.BsmtQual.unique()

array(['Gd', 'TA', 'Ex', 'No Basement', 'Fa'], dtype=object)

dict1={'No Basement':0,'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}

full.BsmtQual=full.BsmtQual.apply(lambda x: dict1[x])

full.BsmtCond.unique()

array(['TA', 'Gd', 'No Basement', 'Fa', 'Po'], dtype=object)

full.BsmtCond=full.BsmtCond.apply(lambda x: dict1[x])

Now let's work with MSZoning.

full.MSZoning.hist()
plt.xlabel("MSZoning")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with Street.

full.Street.hist()
plt.xlabel("Street")
plt.ylabel("Count")
plt.show()

We delete this feature later.

delete.append('Street')

Now let's work with LotShape.

full.LotShape.hist()
plt.xlabel("LotShape")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with LandContour.

full.LandContour.hist()
plt.xlabel("LandContour")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with Utilities.

full.Utilities.hist()
plt.xlabel("Utilities")
plt.ylabel("Count")
plt.show()

full[full.Utilities!="AllPub"]

There is only one entry which has a different value for the feature Utilities.
This feature seems useless.
We will delete this feature later.

delete.append("Utilities")

Now let's work with LotConfig.

full.LotConfig.hist()
plt.xlabel("LotConfig")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with LandSlope.

full.LandSlope.hist()
plt.xlabel("LandSlope")
plt.ylabel("Count")
plt.show()

We will delete this feature later.

delete.append("LandSlope")

Now let's work with Neighborhood.

full.Neighborhood.hist()
plt.xlabel("Neighborhood")
plt.ylabel("Count")
plt.xticks(rotation=90)
plt.show()

We will dummy encode this feature later.

Now let's work with Condition1.

full.Condition1.hist()
plt.xlabel("Condition1")
plt.ylabel("Count")
plt.show()

We will delete this feature later.

delete.append("Condition1")

Now let's work with Condition2.

full.Condition2.hist()
plt.xlabel("Condition2")
plt.ylabel("Count")
plt.show()

We will delete this feature later.

delete.append("Condition2")

Now let's work with BldgType.

full.BldgType.hist()
plt.xlabel("BldgType")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with HouseStyle.

full.HouseStyle.hist()
plt.xlabel("HouseStyle")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with RoofStyle.

full.RoofStyle.hist()
plt.xlabel("RoofStyle")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with RoofMatl.

full.RoofMatl.hist()
plt.xlabel("RoofMatl")
plt.ylabel("Count")
plt.show()

We will delete this feature later.

delete.append("RoofMatl")

Now let's work with Exterior1st.

full.Exterior1st.hist()
plt.xlabel("Exterior1st")
plt.ylabel("Count")
plt.xticks(rotation=90)
plt.show()

We will dummy encode this feature later.

Now let's work with Exterior2nd.

full.Exterior2nd.hist()
plt.xlabel("Exterior2nd")
plt.ylabel("Count")
plt.xticks(rotation=90)
plt.show()

The above two distributions look the same.
We will delete one of them.

delete.append('Exterior2nd')

Now let's work with MasVnrType.

full.MasVnrType.hist()
plt.xlabel("MasVnrType")
plt.ylabel("Count")
plt.xticks(rotation=90)
plt.show()

We will dummy encode this feature later.

Now let's work with ExterQual.

full.ExterQual.hist()
plt.xlabel("ExterQual")
plt.ylabel("Count")
plt.xticks(rotation=90)
plt.show()

dict1={'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}

full.ExterQual=full.ExterQual.apply(lambda x: dict1[x])

Now let's work with ExterCond.

full.ExterCond.hist()
plt.xlabel("ExterCond")
plt.ylabel("Count")
plt.show()

full.ExterCond=full.ExterCond.apply(lambda x: dict1[x])

Now let's work with Foundation.

full.Foundation.hist()
plt.xlabel("Foundation")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with Heating.

full.Heating.hist()
plt.xlabel("Heating")
plt.ylabel("Count")
plt.xticks(rotation=90)
plt.show()

We will remove this feature later.

delete.append("Heating")

Now let's work with HeatingQC.

full.HeatingQC.hist()
plt.xlabel("HeatingQC")
plt.ylabel("Count")
plt.show()

full.HeatingQC=full.HeatingQC.apply(lambda x: dict1[x])

Now let's work with CentralAir.

full.CentralAir.hist()
plt.xlabel("CentralAir")
plt.ylabel("Count")
plt.show()

Let's convert 'Y' to "1" and 'N' to "0". 'Y' is greater than 'N' because with central air you get a higher price.

full.CentralAir=full.CentralAir.apply(lambda x: 1 if x=='Y' else 0)

Now let's work with Electrical.

full.Electrical.hist()
plt.xlabel("Electrical")
plt.ylabel("Count")
plt.show()

We will delete this feature later.

delete.append("Electrical")

Now let's work with KitchenQual.

full.KitchenQual.hist()
plt.xlabel("KitchenQual")
plt.ylabel("Count")
plt.show()

full.KitchenQual=full.KitchenQual.apply(lambda x: dict1[x])

Now let's work with Functional .

full.Functional.hist()
plt.xlabel("Functional")
plt.ylabel("Count")
plt.show()

We will delete this feature later.

delete.append("Functional")

Now let's work with GarageType.

full.GarageType.hist()
plt.xlabel("GarageType")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with PavedDrive.

full.PavedDrive.hist()
plt.xlabel("PavedDrive")
plt.ylabel("Count")
plt.show()

We will delete this feature later.

delete.append("PavedDrive")

Now let's work with SaleType.

full.SaleType.isna().sum()

0

full.SaleType.hist()
plt.xlabel("SaleType")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

Now let's work with SaleCondition.

full.SaleCondition.hist()
plt.xlabel("Sale Condition")
plt.ylabel("Count")
plt.show()

We will dummy encode this feature later.

We are now done with categorical and numerical features.
Let's try adding some new features.

We have some feature that are areas of finished and unfininshed basement and we have a feature that is the total basement array. Let's check if there is a relation.

(full.BsmtFinSF1+full.BsmtFinSF2+full.BsmtUnfSF).equals(full.TotalBsmtSF)

True

We can remove the features BsmtFinSF1, BsmtFinSF2, and BsmtUnfSF.

delete.append("BsmtFinSF1")
delete.append("BsmtFinSF2")
delete.append("BsmtUnfSF")

Let's also check if the sum of 1stFlrSF, 2ndFlrSF, LowQualFinSF is equal to GrLivArea.

(full['1stFlrSF']+full['2ndFlrSF']+full.LowQualFinSF).equals(full.GrLivArea)

True

We can remove the features 1stFlrSF, 2ndFlrSF, and LowQualFinSF.

delete.append("1stFlrSF")
delete.append("2ndFlrSF")
delete.append("LowQualFinSF")

Let's look at the "Bathroom" features now.

print(full.FullBath.corr(y))
print(full.HalfBath.corr(y))
print(full.BsmtFullBath.corr(y))
print(full.BsmtHalfBath.corr(y))

0.5718673954241733
0.2997792902100934
0.2212091220474469
-0.017281473694899883

Let's see if we add the bathroom features together into one feature that may have a better correlation.
Here each half bathroom is multiplied by 0.5.
Credits of this idea goes to Erik Bruin.

full['Total_Bathrooms']=(full.FullBath+full.HalfBath*0.5+full.BsmtFullBath+full.BsmtHalfBath*0.5)

full.Total_Bathrooms.corr(y)

0.6391271749177662

We can remove the other Bath features.

delete.append("FullBath")
delete.append("HalfBath")
delete.append("BsmtFullBath")
delete.append("BsmtHalfBath")

Let's look at the "Porch" features.

print(full.OpenPorchSF.corr(y))
print(full.EnclosedPorch.corr(y))
print(full['3SsnPorch'].corr(y))
print(full.ScreenPorch.corr(y))

0.31554792275247034
-0.11978981907208046
0.02082628480262346
0.11989055290139497

Let's see if adding them improves the correlation with the target feature.

(full.OpenPorchSF+full.EnclosedPorch+full['3SsnPorch']+full.ScreenPorch).corr(y)

0.19782953308159243

It doesn't but I will add them as the correlation in either way is low.

full['TotalPorchSF']=full.OpenPorchSF+full.EnclosedPorch+full['3SsnPorch']+full.ScreenPorch

delete.append("OpenPorchSF")
delete.append("EnclosedPorch")
delete.append("3SsnPorch")
delete.append("ScreenPorch")

full.drop(delete,axis=1,inplace=True)

I will add a feature which describes the age of the house since the last time it was remodeled.

full['Age']=full.YrSold.astype('int')-full.YearRemodAdd.astype('int')

Let's now check the correlation matrix and try to remove highly correlated features.

plt.subplots(figsize=(15,10))
ax=sns.heatmap(full.corr(),cmap=sns.color_palette("RdBu_r", 5),linewidths=.7)
cbar = ax.collections[0].colorbar
cbar.set_ticks(np.arange(-1,1.1,0.2))
plt.title("Heatmap of Numerical Features")
plt.show()

Highly correlated features are chosen whose correlation is shown in red above. We need to delete one of the highly correlated features. To choose one of them, we will compare the correlation of each with the target value and keep the one with the highest correlation.

#compare which of the highly correlated features has a higher correlation with the price
full.PoolArea.corr(y)>full.PoolQC.corr(y)

False

#compare which of the highly correlated features has a higher correlation
full.GarageQual.corr(y)>full.GarageCond.corr(y)

True

#compare which of the highly correlated features has a higher correlation
full.GarageArea.corr(y)>full.GarageCars.corr(y)

False

#compare which of the highly correlated features has a higher correlation
full.ExterQual.corr(y)>full.KitchenQual.corr(y)

True

#drop highly correlated features
full.drop(['PoolArea','GarageCond','GarageArea',"TotRmsAbvGrd",'KitchenQual'],axis=1,inplace=True)

Training

numerical_features,categorical_features=get_features(full)

We are now going to dummy encode the data.
We have an option to drop a column of the dummy variables. This is to ensure we don't fall in the dummy variable trap.

full_wd=full.copy() #wd=with dummy
full_wod=full.copy() #wod= without dummy

full_wd=pd.get_dummies(full_wd)
full_wod=pd.get_dummies(full_wod,drop_first=True)

train_wd=full_wd.iloc[:index,:]
test_wd=full_wd.iloc[index:,:]

train_wod=full_wod.iloc[:index,:]
test_wod=full_wod.iloc[index:,:]

We are going to try several models.
We will install the nedded one here.

!pip install lightgbm

Collecting lightgbm
  Downloading https://files.pythonhosted.org/packages/bf/01/45e209af10fd16537df0c5d8a5474c286554c3eaf9ddb0ce04113f1e8506/lightgbm-2.1.1-py2.py3-none-manylinux1_x86_64.whl (711kB)
    100% |████████████████████████████████| 716kB 1.2MB/s ta 0:00:01
Requirement already satisfied: scipy in /opt/conda/lib/python3.6/site-packages (from lightgbm)
Requirement already satisfied: numpy in /opt/conda/lib/python3.6/site-packages (from lightgbm)
Requirement already satisfied: scikit-learn in /opt/conda/lib/python3.6/site-packages (from lightgbm)
Installing collected packages: lightgbm
Successfully installed lightgbm-2.1.1
You are using pip version 9.0.1, however version 10.0.1 is available.
You should consider upgrading via the 'pip install --upgrade pip' command.

We define the models and assess them here.

Here is our function for calculating the mean rsme and standard deviation of the rsme.

def model_score(name,model,mtrain,y,cv=50):
    
    #mse are by default negative in sklearn
    mse=-cross_val_score(model,mtrain,y,cv=cv,scoring='neg_mean_squared_error')
    rmse=np.sqrt(mse)
    average_rmse=np.mean(rmse)
    std=np.std(rmse)
    #print(rmse)
    print("%s: The average RMSE is: %f with a standard deviation of %f." %(name,average_rmse,std))
    return average_rmse,std

Out of the models that we will use, we have some which are meta-estimators which might give us better results. For more info,click here.

#define models
from sklearn.linear_model import LinearRegression ,ElasticNet, Lasso, LassoLars,HuberRegressor
from sklearn.ensemble import RandomForestRegressor
from sklearn.ensemble import GradientBoostingRegressor,AdaBoostRegressor,BaggingRegressor,ExtraTreesRegressor
from sklearn.svm import  SVR
from lightgbm import LGBMRegressor
from xgboost import XGBRegressor

models=[LinearRegression(),ElasticNet(),Lasso(),LassoLars(),HuberRegressor(),RandomForestRegressor(),\
        GradientBoostingRegressor(),AdaBoostRegressor(),BaggingRegressor(), ExtraTreesRegressor(),SVR(),\
        LGBMRegressor(),XGBRegressor()]


#LR:Linear Regressor
#EN:Elastic Net
#L:Lasso
#Lar:LassoLars
#HR:HuberRegressor
#RF:Random Forest
#GBR:Gradient Boosting Regressor
#Ada:AdaBoostRegressor
#Bag:BaggingRegressor
#Extra:ExtraTreesRegressor
#SVR:Support Vector Regressor
#LGBMR:Light GBM Regressor
#XGBR:XGB Regressor

names=["LR","EN","L","Lar","HR","RF","GBR","Ada","Bag","Extra", "SVR","LGBMR","XGBR"]

Although we don't have to apply k-fold cross validation for ensemble methods such as Random Forests and others since by definition random forests are several decision trees, we will do k-fold cross validation on these models to study the standard deviation of these models.

#lists to save rsme and std of each model
errors=[]
stds=[]
cv_=50
#iterate through models,calculate rsme and std,store them in lists
for name,model in zip(names,models):
    rsme,std=model_score(name,model,train_wod,log_y,cv_)
    errors.append(rsme)
    stds.append(std)
#bar plot of model rsme
bar1=plt.bar(names,errors,label="RSME")
bar2=plt.bar(names,stds,label="Variance")
#we change the color of the bar with min rsme and min std to lighter colors
bar1[np.argmin(errors)].set_color('#bae8f5')
bar2[np.argmin(stds)].set_color('#ffcc66')
plt.legend()

LR: The average RMSE is: 0.128947 with a standard deviation of 0.042698.
EN: The average RMSE is: 0.172177 with a standard deviation of 0.051970.
L: The average RMSE is: 0.179298 with a standard deviation of 0.052386.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.170894 with a standard deviation of 0.052063.
RF: The average RMSE is: 0.148484 with a standard deviation of 0.036599.
GBR: The average RMSE is: 0.127117 with a standard deviation of 0.039404.
Ada: The average RMSE is: 0.174108 with a standard deviation of 0.033708.
Bag: The average RMSE is: 0.148079 with a standard deviation of 0.041187.
Extra: The average RMSE is: 0.146243 with a standard deviation of 0.035813.
SVR: The average RMSE is: 0.396843 with a standard deviation of 0.065563.
LGBMR: The average RMSE is: 0.128632 with a standard deviation of 0.035147.
XGBR: The average RMSE is: 0.126893 with a standard deviation of 0.038967.

<matplotlib.legend.Legend at 0x7f8d04a53198>

#lists to save rsme and std of each model
errors=[]
stds=[]

#iterate through models,calculate rsme and std,store them in lists
for name,model in zip(names,models):
    rsme,std=model_score(name,model,train_wd,log_y)
    errors.append(rsme)
    stds.append(std)
#bar plot of model rsme
bar1=plt.bar(names,errors,label="RSME")
bar2=plt.bar(names,stds,label="Variance")
#we change the color of the bar with min rsme and min std to lighter colors
bar1[np.argmin(errors)].set_color('#bae8f5')
bar2[np.argmin(stds)].set_color('#ffcc66')
plt.legend()

LR: The average RMSE is: 0.128944 with a standard deviation of 0.042698.
EN: The average RMSE is: 0.172177 with a standard deviation of 0.051970.
L: The average RMSE is: 0.179298 with a standard deviation of 0.052386.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.172001 with a standard deviation of 0.052760.
RF: The average RMSE is: 0.145920 with a standard deviation of 0.039825.
GBR: The average RMSE is: 0.127319 with a standard deviation of 0.039875.
Ada: The average RMSE is: 0.174115 with a standard deviation of 0.032819.
Bag: The average RMSE is: 0.144649 with a standard deviation of 0.037005.
Extra: The average RMSE is: 0.143764 with a standard deviation of 0.034898.
SVR: The average RMSE is: 0.396795 with a standard deviation of 0.065595.
LGBMR: The average RMSE is: 0.128511 with a standard deviation of 0.034714.
XGBR: The average RMSE is: 0.126725 with a standard deviation of 0.036599.

<matplotlib.legend.Legend at 0x7f8d048a8978>

Results don't seem better. Feature engineering done was not the best.
Models vary a lot in their results.
In the next approach (and hopefully the final), we will try different version of feature engineering,try to optimize the models by parameter tuning and choose the best one.

Approach #4

#import cleaned data
full=full_cleaned.copy()

Having our cleaned data, we need to select the appropriate features. It turns out that there are some common ways to do so:

RFE: Recursive Forward Elimination
Fitting a random forest and checking the importance of each feature
Fitting a regularized linear model and checking the importance of each feature
Using VIF to eliminate unimportant features.
Trying to get a sense of the data and choosing the features manually.

In this first phase, we will just fit a random forest and a regularized linear model and check the importances of the feature. RFE is a good way to go but we have many features and models to choose from,so we will first try to find the best models and then use RFE to see if we get better results.

Before doing so, we have some features which are considered as numerical but in fact should be categorical.
Let's fix that.

numerical_features,categorical_features=get_features(full)

numerical_features

Index(['MSSubClass', 'LotFrontage', 'LotArea', 'OverallQual', 'OverallCond',
       'YearBuilt', 'YearRemodAdd', 'MasVnrArea', 'BsmtFinSF1', 'BsmtFinSF2',
       'BsmtUnfSF', 'TotalBsmtSF', '1stFlrSF', '2ndFlrSF', 'LowQualFinSF',
       'GrLivArea', 'BsmtFullBath', 'BsmtHalfBath', 'FullBath', 'HalfBath',
       'BedroomAbvGr', 'KitchenAbvGr', 'TotRmsAbvGrd', 'Fireplaces',
       'GarageCars', 'GarageArea', 'WoodDeckSF', 'OpenPorchSF',
       'EnclosedPorch', '3SsnPorch', 'ScreenPorch', 'PoolArea', 'MiscVal',
       'MoSold', 'YrSold'],
      dtype='object')

#MSSubClass should be categorical
full['MSSubClass'] = full['MSSubClass'].apply(str)

#Year and month sold should also be categorical.
full['YrSold'] = full['YrSold'].astype(str)
full['MoSold'] = full['MoSold'].astype(str)

numerical_features,categorical_features=get_features(full)

numerical_features

Index(['LotFrontage', 'LotArea', 'OverallQual', 'OverallCond', 'YearBuilt',
       'YearRemodAdd', 'MasVnrArea', 'BsmtFinSF1', 'BsmtFinSF2', 'BsmtUnfSF',
       'TotalBsmtSF', '1stFlrSF', '2ndFlrSF', 'LowQualFinSF', 'GrLivArea',
       'BsmtFullBath', 'BsmtHalfBath', 'FullBath', 'HalfBath', 'BedroomAbvGr',
       'KitchenAbvGr', 'TotRmsAbvGrd', 'Fireplaces', 'GarageCars',
       'GarageArea', 'WoodDeckSF', 'OpenPorchSF', 'EnclosedPorch', '3SsnPorch',
       'ScreenPorch', 'PoolArea', 'MiscVal'],
      dtype='object')

categorical_features

Index(['MSSubClass', 'MSZoning', 'Street', 'Alley', 'LotShape', 'LandContour',
       'Utilities', 'LotConfig', 'LandSlope', 'Neighborhood', 'Condition1',
       'Condition2', 'BldgType', 'HouseStyle', 'RoofStyle', 'RoofMatl',
       'Exterior1st', 'Exterior2nd', 'MasVnrType', 'ExterQual', 'ExterCond',
       'Foundation', 'BsmtQual', 'BsmtCond', 'BsmtExposure', 'BsmtFinType1',
       'BsmtFinType2', 'Heating', 'HeatingQC', 'CentralAir', 'Electrical',
       'KitchenQual', 'Functional', 'FireplaceQu', 'GarageType',
       'GarageFinish', 'GarageQual', 'GarageCond', 'PavedDrive', 'PoolQC',
       'Fence', 'MiscFeature', 'MoSold', 'YrSold', 'SaleType',
       'SaleCondition'],
      dtype='object')

Let's us encode some categorical features into numbers.
The code is copied form the previous approach and we added a few things.

dict1={'No Pool':0,'Fa':1,'Gd':2,'Ex':3}
full.PoolQC=full.PoolQC.apply(lambda x: dict1[x])

dict1={'NA':0,'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}
full.FireplaceQu=full.FireplaceQu.apply(lambda x: dict1[x])

dict1={'No Garage':0,'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}
full.GarageCond=full.GarageCond.apply(lambda x: dict1[x])
full.GarageQual=full.GarageQual.apply(lambda x: dict1[x])

dict1={'No Garage':0,'Unf':1,'RFn':2,'Fin':3}
full.GarageFinish=full.GarageFinish.apply(lambda x: dict1[x])

dict1={'No Basement':0, 'Unf':1, 'LwQ':2, 'Rec':3, 'BLQ':4, 'ALQ':5, 'GLQ':6}
full.BsmtFinType1=full.BsmtFinType1.apply(lambda x: dict1[x])
full.BsmtFinType2=full.BsmtFinType2.apply(lambda x: dict1[x])

dict1={'No Basement':0, 'No':1, 'Mn':2, 'Av':3, 'Gd':4}
full.BsmtExposure=full.BsmtExposure.apply(lambda x: dict1[x])

dict1={'No Basement':0,'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}
full.BsmtQual=full.BsmtQual.apply(lambda x: dict1[x])
full.BsmtCond=full.BsmtCond.apply(lambda x: dict1[x])

dict1={'Po':1,'TA':2,'Fa':3,'Gd':4,'Ex':5}
full.ExterQual=full.ExterQual.apply(lambda x: dict1[x])
full.ExterCond=full.ExterCond.apply(lambda x: dict1[x])
full.HeatingQC=full.HeatingQC.apply(lambda x: dict1[x])
full.KitchenQual=full.KitchenQual.apply(lambda x: dict1[x])


full.CentralAir=full.CentralAir.apply(lambda x: 1 if x=='Y' else 0)


#added these

#When home functionality is typical,price is higher and it decreases as functionality decreases.
dict1={'Sal':0, 'Sev':1, 'Maj2':2, 'Maj1':3, 'Mod':4, 'Min2':5, 'Min1':6, 'Typ':7}
full.Functional=full.Functional.apply(lambda x: dict1[x])

#Lands with gentle slope usually cost more than lands with moderate slope which usually cost more than lands with severe slopes.
dict1={'Sev':0,'Mod':1,'Gtl':2}
full.LandSlope=full.LandSlope.apply(lambda x: dict1[x])



#A paved street costs more, so let's encode 'Pave' with 1 and 'Grvl' with 0.
full.Street=full.Street.apply(lambda x: 1 if x=='Pave' else 0)

#A paved drive costs more.
dict1={'N':0,'P':1,'Y':2}
full.PavedDrive=full.PavedDrive.apply(lambda x: dict1[x])

Let us also add the features that we added in the last approach.
We will not drop features here,instead we will keep them and will see its importance as a feature for the model later.
We will also add some new features.(credits)

full['Total_Bathrooms']=(full.FullBath+full.HalfBath*0.5+full.BsmtFullBath+full.BsmtHalfBath*0.5)

full['TotalPorchSF']=full.OpenPorchSF+full.EnclosedPorch+full['3SsnPorch']+full.ScreenPorch

#2010 because this dataset was available in 2010
full['Age']=2010-full.YearBuilt

#added these 
full["Remodeled"] = (full["YearRemodAdd"] != full["YearBuilt"]) * 1

area_cols = ['LotFrontage', 'LotArea', 'MasVnrArea', 'BsmtFinSF1', 'BsmtFinSF2', 'BsmtUnfSF',
             'TotalBsmtSF', '1stFlrSF', '2ndFlrSF', 'GrLivArea', 'GarageArea', 'WoodDeckSF', 
             'OpenPorchSF', 'EnclosedPorch', '3SsnPorch', 'ScreenPorch', 'LowQualFinSF', 'PoolArea' ]

full["TotalArea"] =full[area_cols].sum(axis=1)

numerical_features,categorical_features=get_features(full)

But first before fitting the data, we need to encode the categorical variables.

full_hot=pd.get_dummies(full.copy(),drop_first=True)

Now let us see the importance of the features. We will compare feature importances obtained from both lasso model and random forest.

#lasso feature importances
lasso=Lasso()

lasso.fit(full_hot.copy().iloc[:index,:],log_y)
FI_lasso = pd.DataFrame({"Feature Importance":lasso.coef_}, index=full_hot.columns)
FI_lasso[FI_lasso["Feature Importance"]!=0].sort_values("Feature Importance").plot(kind="barh",figsize=(15,25))
plt.show()

#try with random forest
rf=RandomForestRegressor()

rf.fit(full_hot.copy().iloc[:index,:],log_y)
FI_rf = pd.DataFrame({"Feature Importance":rf.feature_importances_}, index=full_hot.columns)
FI_rf.sort_values("Feature Importance").iloc[-20:,:].plot(kind="barh",figsize=(15,25))
plt.show()

The results obtained by Lasso seem weird.
A lot of features that should be there such as OverallQual are not.
On the other hand,the features obtained with the random forest are much more reasonable.
Thus we will take the features obtained with random forest.
The problem here is to choose a certain threshold by which a feature is important.
Inspecting the feature importances values, we choose a threshold of 0.0001. Although the threshold might seem small,but there are much less values and choosing higher thresholds will remove a lot of features.

threshold=0.0001

imp_features=FI_rf[FI_rf['Feature Importance']>threshold].index

full_imp=full_hot.copy()[imp_features]

We now have the most important features.
Let's also extract features using VIF.
The variance inflation factor is a measure for the increase of the variance of the parameter estimates if an additional variable, given by exog_idx is added to the linear regression. It is a measure for multicollinearity of the design matrix, exog.
One recommendation is that if VIF is greater than 5, then the explanatory variable given by exog_idx is highly collinear with the other explanatory variables, and the parameter estimates will have large standard errors because of this.( source)

!pip install statsmodels

Requirement already satisfied: statsmodels in /opt/conda/lib/python3.6/site-packages
You are using pip version 9.0.1, however version 10.0.1 is available.
You should consider upgrading via the 'pip install --upgrade pip' command.

from statsmodels.stats.outliers_influence import variance_inflation_factor

# For each X, calculate VIF and save in dataframe
vif = pd.DataFrame()
vif["VIF Factor"] = [variance_inflation_factor(full_hot.values, i) for i in range(full_hot.shape[1])]
vif["features"] = full_hot.columns

#choose important features
imp_features=vif[vif["VIF Factor"]>5].features.values

full_vif=full_hot.copy()[imp_features]

We have now two different data sets.
We are going to generate more data sets by scaling the data and applying a Box Cox Transformation for the skewed features.
We could have just used the function "log" but applying the box cox transformation gave slightly better results. Credits
For more info about the box cox transformation, click here.

First let's transform the skewed data.
We will use a robust scaler which better deals with outliers. From these 2 datasets, we will produce another 2 where we apply the box cox transformation for skewed variables.

from scipy.special import boxcox1p

#for full_imp
skew=pd.DataFrame()
skew['Skew']=full_imp.skew()
skewed_columns=skew[skew>0.7].dropna().index
full_imp_skewed=full_imp.copy()
full_imp_skewed[skewed_columns]=boxcox1p(full_imp_skewed[skewed_columns],0.15)


#for full_vif
skew=pd.DataFrame()
skew['Skew']=full_vif.skew()
skewed_columns=skew[skew>0.7].dropna().index
full_vif_skewed=full_vif.copy()
full_vif_skewed[skewed_columns]=boxcox1p(full_vif_skewed[skewed_columns],0.15)

Now let's transform the skewed data.
We will use a robust scaler which better deals with outliers. From the 4 datasets we have, we will produce another 4 where are a scaled version.

from sklearn.preprocessing import RobustScaler

scaler=RobustScaler()

full_imp_scaled=full_imp.copy()
full_imp_scaled[full_imp.columns]=scaler.fit_transform(full_imp_scaled)

full_imp_skewed_scaled=full_imp_skewed.copy()
full_imp_skewed_scaled[full_imp_skewed.columns]=scaler.fit_transform(full_imp_skewed_scaled)

full_vif_scaled=full_vif.copy()
full_vif_scaled[full_vif_scaled.columns]=scaler.fit_transform(full_vif_scaled)

full_vif_skewed_scaled=full_vif_skewed.copy()
full_vif_skewed_scaled[full_vif_skewed_scaled.columns]=scaler.fit_transform(full_vif_skewed_scaled)

We now have 8 different data sets.
Let's join them in one list.

datasets=[full_imp,full_imp_skewed,full_imp_scaled,full_imp_skewed_scaled,full_vif,full_vif_skewed,full_vif_scaled,full_vif_skewed_scaled]

datasets_names=["full_imp","full_imp_skewed","full_imp_scaled","full_imp_skewed_scaled","full_vif","full_vif_skewed","full_vif_scaled","full_vif_skewed_scaled"]

We can now train our models.

from sklearn.cross_validation import cross_val_score
def model_score(name,model,mtrain,y,cv=50):
    
    #mse are by default negative in sklearn
    mse=-cross_val_score(model,mtrain,y,cv=cv,scoring='neg_mean_squared_error')
    rmse=np.sqrt(mse)
    average_rmse=np.mean(rmse)
    std=np.std(rmse)
    #print(rmse)
    print("%s: The average RMSE is: %f with a standard deviation of %f." %(name,average_rmse,std))
    return average_rmse,std

#define models
from sklearn.linear_model import LinearRegression ,ElasticNet, Lasso, LassoLars,HuberRegressor,Ridge,BayesianRidge
from sklearn.kernel_ridge import KernelRidge
from sklearn.ensemble import RandomForestRegressor
from sklearn.ensemble import GradientBoostingRegressor,AdaBoostRegressor,BaggingRegressor,ExtraTreesRegressor
from sklearn.svm import  SVR
from lightgbm import LGBMRegressor
from xgboost import XGBRegressor

models=[LinearRegression(),ElasticNet(),Lasso(),LassoLars(),HuberRegressor(),RandomForestRegressor(),\
        GradientBoostingRegressor(n_estimators=500,loss='huber'),AdaBoostRegressor(),BaggingRegressor(), ExtraTreesRegressor(),SVR(),\
        LGBMRegressor(objective='regression',n_estimators=500),XGBRegressor(n_estimators=500),KernelRidge(kernel='polynomial'),Ridge(),BayesianRidge()]


#LR:Linear Regressor
#EN:Elastic Net
#L:Lasso
#Lar:LassoLars
#HR:HuberRegressor
#RF:Random Forest
#GBR:Gradient Boosting Regressor
#Ada:AdaBoostRegressor
#Bag:BaggingRegressor
#Extra:ExtraTreesRegressor
#SVR:Support Vector Regressor
#LGBMR:Light GBM Regressor
#XGBR:XGB Regressor
#KR: Kernel Ridge
#R:Ridge
#BR:Bayesian Ridge
names=["LR","EN","L","Lar","HR","RF","GBR","Ada","Bag","Extra", "SVR","LGBMR","XGBR","KR","R","BR"]

for i in range(len(datasets)):
    data=datasets[i]
    
    #get train dataset
    train=data.iloc[:index,:]

    #lists to save rsme and std of each model
    errors=[]
    stds=[]
    cv_=50
    #iterate through models,calculate rsme and std,store them in lists
    for name,model in zip(names,models):
        rsme,std=model_score(name,model,train,log_y,cv_)
        errors.append(rsme)
        stds.append(std)
    #bar plot of model rsme
    plt.figure()
    plt.title(datasets_names[i])
    plt.ylim([0,0.5])
    bar1=plt.bar(names,errors,label="RSME")
    bar2=plt.bar(names,stds,label="Variance")
    #we change the color of the bar with min rsme and min std to lighter colors
    bar1[np.argmin(errors)].set_color('#bae8f5')
    bar2[np.argmin(stds)].set_color('#ffcc66')
    plt.legend()
    plt.show()

LR: The average RMSE is: 0.124118 with a standard deviation of 0.045235.
EN: The average RMSE is: 0.163913 with a standard deviation of 0.051284.
L: The average RMSE is: 0.171157 with a standard deviation of 0.052029.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.162140 with a standard deviation of 0.051853.
RF: The average RMSE is: 0.143477 with a standard deviation of 0.036276.
GBR: The average RMSE is: 0.123962 with a standard deviation of 0.041423.
Ada: The average RMSE is: 0.170529 with a standard deviation of 0.034180.
Bag: The average RMSE is: 0.143416 with a standard deviation of 0.040869.
Extra: The average RMSE is: 0.142532 with a standard deviation of 0.035357.
SVR: The average RMSE is: 0.397461 with a standard deviation of 0.065301.
LGBMR: The average RMSE is: 0.124818 with a standard deviation of 0.033953.
XGBR: The average RMSE is: 0.124171 with a standard deviation of 0.038004.
KR: The average RMSE is: 5.031851 with a standard deviation of 22.653322.
R: The average RMSE is: 0.122880 with a standard deviation of 0.045303.
BR: The average RMSE is: 0.120183 with a standard deviation of 0.045882.

LR: The average RMSE is: 0.122452 with a standard deviation of 0.036513.
EN: The average RMSE is: 0.259958 with a standard deviation of 0.049730.
L: The average RMSE is: 0.263615 with a standard deviation of 0.051129.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.282480 with a standard deviation of 0.065255.
RF: The average RMSE is: 0.148326 with a standard deviation of 0.041908.
GBR: The average RMSE is: 0.124987 with a standard deviation of 0.041626.
Ada: The average RMSE is: 0.170051 with a standard deviation of 0.032370.
Bag: The average RMSE is: 0.145537 with a standard deviation of 0.037147.
Extra: The average RMSE is: 0.142720 with a standard deviation of 0.037727.
SVR: The average RMSE is: 0.319470 with a standard deviation of 0.065922.
LGBMR: The average RMSE is: 0.124807 with a standard deviation of 0.033930.
XGBR: The average RMSE is: 0.124069 with a standard deviation of 0.038015.
KR: The average RMSE is: 0.296542 with a standard deviation of 0.114131.
R: The average RMSE is: 0.120369 with a standard deviation of 0.036691.
BR: The average RMSE is: 0.118133 with a standard deviation of 0.036353.

LR: The average RMSE is: 0.124220 with a standard deviation of 0.045304.
EN: The average RMSE is: 0.392857 with a standard deviation of 0.064027.
L: The average RMSE is: 0.393253 with a standard deviation of 0.063663.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 1.458720 with a standard deviation of 0.278596.
RF: The average RMSE is: 0.147080 with a standard deviation of 0.038503.
GBR: The average RMSE is: 0.124353 with a standard deviation of 0.042702.
Ada: The average RMSE is: 0.170995 with a standard deviation of 0.033179.
Bag: The average RMSE is: 0.147939 with a standard deviation of 0.037449.
Extra: The average RMSE is: 0.144019 with a standard deviation of 0.036886.
SVR: The average RMSE is: 0.195149 with a standard deviation of 0.063452.
LGBMR: The average RMSE is: 0.125113 with a standard deviation of 0.033151.
XGBR: The average RMSE is: 0.124279 with a standard deviation of 0.038021.
KR: The average RMSE is: 1.105860 with a standard deviation of 1.339801.
R: The average RMSE is: 0.122721 with a standard deviation of 0.045285.
BR: The average RMSE is: 0.119941 with a standard deviation of 0.045148.

LR: The average RMSE is: 0.125612 with a standard deviation of 0.039048.
EN: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
L: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.145000 with a standard deviation of 0.042663.
RF: The average RMSE is: 0.146933 with a standard deviation of 0.038445.
GBR: The average RMSE is: 0.123915 with a standard deviation of 0.043270.
Ada: The average RMSE is: 0.173475 with a standard deviation of 0.034072.
Bag: The average RMSE is: 0.144262 with a standard deviation of 0.039097.
Extra: The average RMSE is: 0.141477 with a standard deviation of 0.038251.
SVR: The average RMSE is: 0.128022 with a standard deviation of 0.041803.
LGBMR: The average RMSE is: 0.125096 with a standard deviation of 0.033130.
XGBR: The average RMSE is: 0.124178 with a standard deviation of 0.038033.
KR: The average RMSE is: 0.133004 with a standard deviation of 0.043544.
R: The average RMSE is: 0.120311 with a standard deviation of 0.036709.
BR: The average RMSE is: 0.117944 with a standard deviation of 0.036482.

LR: The average RMSE is: 190.731013 with a standard deviation of 1080.517915.
EN: The average RMSE is: 0.174541 with a standard deviation of 0.053177.
L: The average RMSE is: 0.177836 with a standard deviation of 0.052648.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.175643 with a standard deviation of 0.053243.
RF: The average RMSE is: 0.151557 with a standard deviation of 0.038649.
GBR: The average RMSE is: 0.142244 with a standard deviation of 0.042493.
Ada: The average RMSE is: 0.178623 with a standard deviation of 0.035352.
Bag: The average RMSE is: 0.153763 with a standard deviation of 0.039021.
Extra: The average RMSE is: 0.150416 with a standard deviation of 0.039792.
SVR: The average RMSE is: 0.397433 with a standard deviation of 0.065311.
LGBMR: The average RMSE is: 0.145108 with a standard deviation of 0.033010.
XGBR: The average RMSE is: 0.139838 with a standard deviation of 0.040911.
KR: The average RMSE is: 2.857032 with a standard deviation of 10.161093.
R: The average RMSE is: 0.138667 with a standard deviation of 0.043520.
BR: The average RMSE is: 0.135393 with a standard deviation of 0.045562.

LR: The average RMSE is: 3603.501612 with a standard deviation of 22687.128280.
EN: The average RMSE is: 0.275552 with a standard deviation of 0.052568.
L: The average RMSE is: 0.277220 with a standard deviation of 0.052972.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.261023 with a standard deviation of 0.070999.
RF: The average RMSE is: 0.153671 with a standard deviation of 0.039764.
GBR: The average RMSE is: 0.143295 with a standard deviation of 0.044679.
Ada: The average RMSE is: 0.178823 with a standard deviation of 0.033318.
Bag: The average RMSE is: 0.151137 with a standard deviation of 0.039644.
Extra: The average RMSE is: 0.148088 with a standard deviation of 0.038925.
SVR: The average RMSE is: 0.295764 with a standard deviation of 0.060064.
LGBMR: The average RMSE is: 0.145044 with a standard deviation of 0.032920.
XGBR: The average RMSE is: 0.139805 with a standard deviation of 0.040917.
KR: The average RMSE is: 0.537151 with a standard deviation of 0.224184.
R: The average RMSE is: 0.134610 with a standard deviation of 0.035835.
BR: The average RMSE is: 0.132392 with a standard deviation of 0.036994.

LR: The average RMSE is: 0.145249 with a standard deviation of 0.046849.
EN: The average RMSE is: 0.398127 with a standard deviation of 0.077649.
L: The average RMSE is: 0.397766 with a standard deviation of 0.074542.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 2.150668 with a standard deviation of 0.381973.
RF: The average RMSE is: 0.154684 with a standard deviation of 0.039999.
GBR: The average RMSE is: 0.142801 with a standard deviation of 0.044055.
Ada: The average RMSE is: 0.178569 with a standard deviation of 0.033476.
Bag: The average RMSE is: 0.156491 with a standard deviation of 0.037155.
Extra: The average RMSE is: 0.149674 with a standard deviation of 0.040028.
SVR: The average RMSE is: 0.218709 with a standard deviation of 0.062930.
LGBMR: The average RMSE is: 0.146080 with a standard deviation of 0.032177.
XGBR: The average RMSE is: 0.139831 with a standard deviation of 0.040903.
KR: The average RMSE is: 103.532565 with a standard deviation of 698.953681.
R: The average RMSE is: 0.138670 with a standard deviation of 0.043473.
BR: The average RMSE is: 0.135543 with a standard deviation of 0.044571.

LR: The average RMSE is: 11475979.424061 with a standard deviation of 50230689.945059.
EN: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
L: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
Lar: The average RMSE is: 0.398519 with a standard deviation of 0.064352.
HR: The average RMSE is: 0.203946 with a standard deviation of 0.072745.
RF: The average RMSE is: 0.154108 with a standard deviation of 0.040019.
GBR: The average RMSE is: 0.142570 with a standard deviation of 0.044024.
Ada: The average RMSE is: 0.177743 with a standard deviation of 0.034859.
Bag: The average RMSE is: 0.154938 with a standard deviation of 0.039456.
Extra: The average RMSE is: 0.151775 with a standard deviation of 0.037277.
SVR: The average RMSE is: 0.147555 with a standard deviation of 0.044955.
LGBMR: The average RMSE is: 0.146068 with a standard deviation of 0.032149.
XGBR: The average RMSE is: 0.139805 with a standard deviation of 0.040917.
KR: The average RMSE is: 0.145806 with a standard deviation of 0.044280.
R: The average RMSE is: 0.134545 with a standard deviation of 0.035828.
BR: The average RMSE is: 0.132278 with a standard deviation of 0.036618.

It seems that features using the Random Forest feature importance have produced beter results than using VIF.
The Bayesian Ridge produced the lowest error in all datasets.

The lowest result was using the data set "full_imp_skewed_scaled",so we will from now on train the models in this data set.

Let's try doing some grid search to try and optimize the models.
We will choose the models with the lowest errors in this data set which are Bayesian Ridge,Ridge,Gradient Boosting Regressor, and XGBoost.

from sklearn.model_selection import GridSearchCV

We start with the Bayesian Ridge.

train=full_imp_skewed_scaled.iloc[:index,:]

parameters=[{'n_iter':[300,1000],'tol':[0.1,0.001],'alpha_1':[1e-6,1e-4],\
             'lambda_1':[1e-6,1e-4]}]

BR=BayesianRidge()
grid_search= GridSearchCV(estimator=BR,param_grid=parameters,scoring='neg_mean_squared_error',cv=50,n_jobs=-1,verbose=1)
grid_search=grid_search.fit(train,log_y)
print("The best score is %f." %(grid_search.best_score_))
print("The best parameters are ",grid_search.best_params_ )
grid_search.grid_scores_

Fitting 50 folds for each of 16 candidates, totalling 800 fits

[Parallel(n_jobs=-1)]: Done 136 tasks      | elapsed:   29.4s
[Parallel(n_jobs=-1)]: Done 386 tasks      | elapsed:  1.7min
[Parallel(n_jobs=-1)]: Done 736 tasks      | elapsed:  3.1min

The best score is -0.015141.
The best parameters are  {'alpha_1': 1e-06, 'n_iter': 300, 'lambda_1': 0.0001, 'tol': 0.001}

[Parallel(n_jobs=-1)]: Done 800 out of 800 | elapsed:  3.3min finished
/usr/local/lib/python3.5/dist-packages/sklearn/model_selection/_search.py:761: DeprecationWarning: The grid_scores_ attribute was deprecated in version 0.18 in favor of the more elaborate cv_results_ attribute. The grid_scores_ attribute will not be available from 0.20
  DeprecationWarning)

[mean: -0.01515, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 300, 'lambda_1': 1e-06, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 300, 'lambda_1': 1e-06, 'tol': 0.001},
 mean: -0.01515, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 1000, 'lambda_1': 1e-06, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 1000, 'lambda_1': 1e-06, 'tol': 0.001},
 mean: -0.01515, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 300, 'lambda_1': 0.0001, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 300, 'lambda_1': 0.0001, 'tol': 0.001},
 mean: -0.01515, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 1000, 'lambda_1': 0.0001, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 1e-06, 'n_iter': 1000, 'lambda_1': 0.0001, 'tol': 0.001},
 mean: -0.01515, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 300, 'lambda_1': 1e-06, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 300, 'lambda_1': 1e-06, 'tol': 0.001},
 mean: -0.01515, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 1000, 'lambda_1': 1e-06, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 1000, 'lambda_1': 1e-06, 'tol': 0.001},
 mean: -0.01515, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 300, 'lambda_1': 0.0001, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 300, 'lambda_1': 0.0001, 'tol': 0.001},
 mean: -0.01515, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 1000, 'lambda_1': 0.0001, 'tol': 0.1},
 mean: -0.01514, std: 0.01178, params: {'alpha_1': 0.0001, 'n_iter': 1000, 'lambda_1': 0.0001, 'tol': 0.001}]

b=BayesianRidge(alpha_1=1e-6,n_iter=300,lambda_1=0.0001,tol=0.001)
mse=-cross_val_score(b,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.118915 with a standard deviation of 0.035495.

Let's look at Ridge now.

parameters=[{'alpha':[1.e-4,0.1,1,5,10],'fit_intercept':[True,False],'tol':[0.001,0.1,1]}]

R=Ridge()
grid_search= GridSearchCV(estimator=R,param_grid=parameters,scoring='neg_mean_squared_error',cv=50,n_jobs=-1,verbose=1)
grid_search=grid_search.fit(train,log_y)
print("The best score is %f." %(grid_search.best_score_))
print("The best parameters are ",grid_search.best_params_ )

Fitting 50 folds for each of 30 candidates, totalling 1500 fits

[Parallel(n_jobs=-1)]: Done 199 tasks      | elapsed:    9.4s
[Parallel(n_jobs=-1)]: Done 449 tasks      | elapsed:   25.7s
[Parallel(n_jobs=-1)]: Done 799 tasks      | elapsed:   48.0s
[Parallel(n_jobs=-1)]: Done 1249 tasks      | elapsed:  1.3min

The best score is -0.015119.
The best parameters are  {'alpha': 10, 'fit_intercept': True, 'tol': 0.001}

[Parallel(n_jobs=-1)]: Done 1500 out of 1500 | elapsed:  1.6min finished

R=Ridge(alpha=10,fit_intercept=True,tol=0.001)
mse=-cross_val_score(R,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.118841 with a standard deviation of 0.035453.

Let's look at Gradient Boosting Regressor.

parameters=[{'loss':['ls','huber'],'learning_rate':[0.1,0.01,1],'n_estimators':[100,300,500],'max_depth':[3,7,10,15]}]

GBR=GradientBoostingRegressor()
grid_search= GridSearchCV(estimator=GBR,param_grid=parameters,scoring='neg_mean_squared_error',cv=50,n_jobs=-1,verbose=1)
grid_search=grid_search.fit(train,log_y)
print("The best score is %f." %(grid_search.best_score_))
print("The best parameters are ",grid_search.best_params_ )

Fitting 50 folds for each of 72 candidates, totalling 3600 fits

[Parallel(n_jobs=-1)]: Done 136 tasks      | elapsed:   31.4s
[Parallel(n_jobs=-1)]: Done 386 tasks      | elapsed:  3.5min
[Parallel(n_jobs=-1)]: Done 736 tasks      | elapsed:  6.7min
[Parallel(n_jobs=-1)]: Done 1186 tasks      | elapsed: 15.2min
[Parallel(n_jobs=-1)]: Done 1736 tasks      | elapsed: 31.4min
[Parallel(n_jobs=-1)]: Done 2386 tasks      | elapsed: 69.4min
[Parallel(n_jobs=-1)]: Done 3136 tasks      | elapsed: 72.2min
[Parallel(n_jobs=-1)]: Done 3600 out of 3600 | elapsed: 73.6min finished

The best score is -0.016856.
The best parameters are  {'n_estimators': 500, 'learning_rate': 0.1, 'max_depth': 3, 'loss': 'ls'}

GBR=GradientBoostingRegressor(loss='ls',n_estimators=500,learning_rate=0.1,max_depth=3)
mse=-cross_val_score(GBR,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.124749 with a standard deviation of 0.040219.

Let's look at XGBoost.

parameters=[{'booster':['gbtree','dart'],'n_estimators':[100,300],'gamma':[0,0.5,1],'max_depth':[3,7,10],\
             'reg_lambda':[0.1,1],'reg_alpha':[0.1,1]}]

XGBR=XGBRegressor()
grid_search= GridSearchCV(estimator=XGBR,param_grid=parameters,scoring='neg_mean_squared_error',cv=50,n_jobs=-1,verbose=1)
grid_search=grid_search.fit(train,log_y)
print("The best score is %f." %(grid_search.best_score_))
print("The best parameters are ",grid_search.best_params_ )
grid_search.grid_scores_

Fitting 50 folds for each of 144 candidates, totalling 7200 fits

[Parallel(n_jobs=-1)]: Done 136 tasks      | elapsed:   16.9s
[Parallel(n_jobs=-1)]: Done 386 tasks      | elapsed:  1.3min
[Parallel(n_jobs=-1)]: Done 736 tasks      | elapsed:  3.2min
[Parallel(n_jobs=-1)]: Done 1186 tasks      | elapsed:  6.2min
[Parallel(n_jobs=-1)]: Done 1736 tasks      | elapsed:  8.7min
[Parallel(n_jobs=-1)]: Done 2386 tasks      | elapsed: 15.4min
[Parallel(n_jobs=-1)]: Done 3136 tasks      | elapsed: 19.1min
[Parallel(n_jobs=-1)]: Done 3986 tasks      | elapsed: 27.4min
[Parallel(n_jobs=-1)]: Done 4936 tasks      | elapsed: 43.8min
[Parallel(n_jobs=-1)]: Done 5986 tasks      | elapsed: 58.5min
[Parallel(n_jobs=-1)]: Done 7136 tasks      | elapsed: 73.3min
[Parallel(n_jobs=-1)]: Done 7200 out of 7200 | elapsed: 74.6min finished

The best score is -0.016699.
The best parameters are  {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0}

/usr/local/lib/python3.5/dist-packages/sklearn/model_selection/_search.py:761: DeprecationWarning: The grid_scores_ attribute was deprecated in version 0.18 in favor of the more elaborate cv_results_ attribute. The grid_scores_ attribute will not be available from 0.20
  DeprecationWarning)

[mean: -0.01787, std: 0.01214, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01793, std: 0.01138, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01822, std: 0.01056, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01862, std: 0.01114, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01720, std: 0.01209, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01704, std: 0.01135, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01671, std: 0.01050, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01695, std: 0.01072, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01830, std: 0.01195, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01818, std: 0.01135, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01820, std: 0.01003, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01886, std: 0.01067, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01825, std: 0.01202, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01799, std: 0.01128, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01782, std: 0.00992, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01855, std: 0.01053, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01866, std: 0.01187, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01819, std: 0.01101, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01909, std: 0.01050, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01937, std: 0.01104, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01862, std: 0.01191, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01817, std: 0.01107, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01880, std: 0.01044, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01907, std: 0.01090, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.02142, std: 0.01132, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02192, std: 0.01142, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02228, std: 0.01193, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02222, std: 0.01168, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02142, std: 0.01132, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02192, std: 0.01142, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02228, std: 0.01193, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02221, std: 0.01167, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02133, std: 0.01131, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02133, std: 0.01131, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02126, std: 0.01123, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02126, std: 0.01124, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02400, std: 0.01240, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02415, std: 0.01204, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02473, std: 0.01266, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02415, std: 0.01204, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02473, std: 0.01266, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'gbtree', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'gbtree', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.01782, std: 0.01214, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01793, std: 0.01138, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01822, std: 0.01056, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01862, std: 0.01114, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01711, std: 0.01205, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01703, std: 0.01122, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01670, std: 0.01052, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01697, std: 0.01070, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01824, std: 0.01195, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01815, std: 0.01132, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01819, std: 0.01002, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01888, std: 0.01066, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01812, std: 0.01204, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01794, std: 0.01127, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01781, std: 0.00995, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01859, std: 0.01052, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01868, std: 0.01188, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01822, std: 0.01109, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01903, std: 0.01041, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01939, std: 0.01104, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0},
 mean: -0.01863, std: 0.01189, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01819, std: 0.01114, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01876, std: 0.01040, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.01907, std: 0.01091, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0},
 mean: -0.02142, std: 0.01132, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02192, std: 0.01142, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02228, std: 0.01193, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02222, std: 0.01168, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02142, std: 0.01132, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02192, std: 0.01142, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02228, std: 0.01193, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02221, std: 0.01167, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02133, std: 0.01131, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02133, std: 0.01131, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02126, std: 0.01123, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 0.5},
 mean: -0.02126, std: 0.01124, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02183, std: 0.01156, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02232, std: 0.01201, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02225, std: 0.01168, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 0.5},
 mean: -0.02400, std: 0.01240, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02415, std: 0.01204, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02473, std: 0.01266, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02415, std: 0.01204, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02473, std: 0.01266, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 3, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 7, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 100, 'gamma': 1},
 mean: -0.02400, std: 0.01240, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02411, std: 0.01212, params: {'booster': 'dart', 'reg_alpha': 0.1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02475, std: 0.01265, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 0.1, 'n_estimators': 300, 'gamma': 1},
 mean: -0.02485, std: 0.01289, params: {'booster': 'dart', 'reg_alpha': 1, 'max_depth': 10, 'reg_lambda': 1, 'n_estimators': 300, 'gamma': 1}]

XGBR=XGBRegressor(booster='dart',reg_alpha=1,max_depth=3,reg_lambda=0.1,n_estimators=300,gamma=0)
mse=-cross_val_score(XGBR,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.125201 with a standard deviation of 0.035950.

We see that we get slightly better results/
Let us now try building a better regressor by stacking several regressors together.
We write a simple class that takes some models and returns the average of their predicitons.

from sklearn.base import BaseEstimator, TransformerMixin, RegressorMixin, clone

#https://www.kaggle.com/serigne/stacked-regressions-top-4-on-leaderboard

class AveragingModels(BaseEstimator, RegressorMixin, TransformerMixin):
    def __init__(self, models):
        self.models = models
        
    # we define clones of the original models to fit the data in
    def fit(self, X, y):
        self.models_ = [clone(x) for x in self.models]
        
        # Train cloned base models
        for model in self.models_:
            model.fit(X, y)

        return self
    
    #Now we do the predictions for cloned models and average them
    def predict(self, X):
        predictions = np.column_stack([
            model.predict(X) for model in self.models_
        ])
        return np.mean(predictions, axis=1)

#average
averaged_models = AveragingModels(models = (b, R, GBR, XGBR))

mse=-cross_val_score(averaged_models,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.113766 with a standard deviation of 0.036917.

Now I am going to weight the models arbitrarly and see if we get better results.

#https://www.kaggle.com/serigne/stacked-regressions-top-4-on-leaderboard

class AveragingModels(BaseEstimator, RegressorMixin, TransformerMixin):
    def __init__(self, models):
        self.models = models
        
    # we define clones of the original models to fit the data in
    def fit(self, X, y):
        self.models_ = [clone(x) for x in self.models]
        
        # Train cloned base models
        for model in self.models_:
            model.fit(X, y)

        return self
    
    #Now we do the predictions for cloned models and average them
    def predict(self, X):
        predictions = np.column_stack([
            model.predict(X) for model in self.models_
        ])
        
        w=np.array([0.3,0.3,0.2,0.2])
        
        return np.dot(predictions, w.T)

#average
averaged_models = AveragingModels(models = (b, R, GBR, XGBR))

mse=-cross_val_score(averaged_models,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.113737 with a standard deviation of 0.036959.

let's try other values for the weights.

#https://www.kaggle.com/serigne/stacked-regressions-top-4-on-leaderboard

class AveragingModels(BaseEstimator, RegressorMixin, TransformerMixin):
    def __init__(self, models):
        self.models = models
        
    # we define clones of the original models to fit the data in
    def fit(self, X, y):
        self.models_ = [clone(x) for x in self.models]
        
        # Train cloned base models
        for model in self.models_:
            model.fit(X, y)

        return self
    
    #Now we do the predictions for cloned models and average them
    def predict(self, X):
        predictions = np.column_stack([
            model.predict(X) for model in self.models_
        ])
        
        w=np.array([0.4,0.3,0.2,0.1])
        
        return np.dot(predictions, w.T)

#average
averaged_models = AveragingModels(models = (b, R, GBR, XGBR))

mse=-cross_val_score(averaged_models,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.113836 with a standard deviation of 0.036795.

Let's try one more.

# credits: https://www.kaggle.com/serigne/stacked-regressions-top-4-on-leaderboard

class AveragingModels(BaseEstimator, RegressorMixin, TransformerMixin):
    def __init__(self, models):
        self.models = models
        
    # we define clones of the original models to fit the data in
    def fit(self, X, y):
        self.models_ = [clone(x) for x in self.models]
        
        # Train cloned base models
        for model in self.models_:
            model.fit(X, y)

        return self
    
    #Now we do the predictions for cloned models and average them
    def predict(self, X):
        predictions = np.column_stack([
            model.predict(X) for model in self.models_
        ])
        
        w=np.array([0.5,0.2,0.1,0.1])
        
        return np.dot(predictions, w.T)

#average
averaged_models = AveragingModels(models = (b, R, GBR, XGBR))

mse=-cross_val_score(averaged_models,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 1.209081 with a standard deviation of 0.024228.

Let's try a Neural Network using Keras.

!pip install keras

Collecting keras
  Downloading https://files.pythonhosted.org/packages/54/e8/eaff7a09349ae9bd40d3ebaf028b49f5e2392c771f294910f75bb608b241/Keras-2.1.6-py2.py3-none-any.whl (339kB)
    100% |################################| 348kB 1.6MB/s ta 0:00:01
Requirement already satisfied: h5py in /usr/local/lib/python3.5/dist-packages (from keras)
Collecting pyyaml (from keras)
  Downloading https://files.pythonhosted.org/packages/4a/85/db5a2df477072b2902b0eb892feb37d88ac635d36245a72a6a69b23b383a/PyYAML-3.12.tar.gz (253kB)
    100% |################################| 256kB 2.0MB/s ta 0:00:01
Requirement already satisfied: numpy>=1.9.1 in /usr/local/lib/python3.5/dist-packages (from keras)
Requirement already satisfied: scipy>=0.14 in /usr/local/lib/python3.5/dist-packages (from keras)
Requirement already satisfied: six>=1.9.0 in /usr/local/lib/python3.5/dist-packages (from keras)
Building wheels for collected packages: pyyaml
  Running setup.py bdist_wheel for pyyaml ... done
  Stored in directory: /root/.cache/pip/wheels/03/05/65/bdc14f2c6e09e82ae3e0f13d021e1b6b2481437ea2f207df3f
Successfully built pyyaml
Installing collected packages: pyyaml, keras
Successfully installed keras-2.1.6 pyyaml-3.12
You are using pip version 9.0.1, however version 10.0.1 is available.
You should consider upgrading via the 'pip install --upgrade pip' command.

from keras.models import Sequential
from keras.layers import Dense
from keras.wrappers.scikit_learn import KerasRegressor
from sklearn.model_selection import cross_val_score
from sklearn.model_selection import KFold
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import Pipeline

def baseline_model():
    model = Sequential()
    model.add(Dense(130, input_dim=130,kernel_initializer='normal', activation='relu'))
    model.add(Dense(80, activation='relu',kernel_initializer='normal'))
    model.add(Dense(40, activation='relu',kernel_initializer='normal'))
    model.add(Dense(20, activation='relu',kernel_initializer='normal'))
    model.add(Dense(10, activation='relu',kernel_initializer='normal'))
    model.add(Dense(1,kernel_initializer='normal'))
    # Compile model
    model.compile(loss='mean_squared_error', optimizer='adam')
    return model

# fix random seed for reproducibility
seed = 7
numpy.random.seed(seed)

# evaluate model
estimator = KerasRegressor(build_fn=baseline_model, epochs=100, batch_size=5, verbose=0)

kfold = KFold(n_splits=50, random_state=seed)
results = cross_val_score(estimator, train, log_y, cv=kfold,verbose=1)
print("Results: %.2f (%.2f) MSE" % (results.mean(), results.std()))

Results: -0.17 (0.09) MSE

[Parallel(n_jobs=1)]: Done  50 out of  50 | elapsed: 54.5min finished

The best model seems to be the stacked one.

class AveragingModels(BaseEstimator, RegressorMixin, TransformerMixin):
    def __init__(self, models):
        self.models = models
        
    # we define clones of the original models to fit the data in
    def fit(self, X, y):
        self.models_ = [clone(x) for x in self.models]
        
        # Train cloned base models
        for model in self.models_:
            model.fit(X, y)

        return self
    
    #Now we do the predictions for cloned models and average them
    def predict(self, X):
        #self.models_ = [clone(x) for x in self.models]
        predictions = np.column_stack([
            model.predict(X) for model in self.models_
        ])
        
        w=np.array([0.3,0.3,0.2,0.2])
        
        return np.dot(predictions, w.T)

#average
averaged_models = AveragingModels(models = (b, R, GBR, XGBR))
mse=-cross_val_score(averaged_models,train,log_y,cv=50,scoring='neg_mean_squared_error')
rmse=np.sqrt(mse)
average_rmse=np.mean(rmse)
std=np.std(rmse)
#print(rmse)
print(" The average RMSE is: %f with a standard deviation of %f." %(average_rmse,std))

 The average RMSE is: 0.113683 with a standard deviation of 0.037014.

test=full_imp_skewed_scaled.iloc[index:,:]

test.shape

(260, 130)

averaged_models.fit(train,log_y)
y_pred=averaged_models.predict(test)

y_pred=np.exp(y_pred)

submissions=pd.DataFrame({'Id':test_ID,'SalePrice':y_pred})
submissions.to_csv('submission.csv',index=False)

Before we conclude,we want to test another neural network. The following neural network is just for comparison. Credits for the code go entirely to Johnny Liu. He uses a package called Gluon from the MXNet library.

!pip install mxnet

Collecting mxnet
  Downloading https://files.pythonhosted.org/packages/96/98/c9877e100c3d1ac92263bfaba7bb8a49294e099046592040a2ff8620ac61/mxnet-1.1.0.post0-py2.py3-none-manylinux1_x86_64.whl (23.8MB)
    100% |████████████████████████████████| 23.8MB 58kB/s  eta 0:00:01
Requirement already satisfied: numpy<1.15.0,>=1.8.2 in /opt/conda/lib/python3.6/site-packages (from mxnet)
Collecting graphviz<0.9.0,>=0.8.1 (from mxnet)
  Downloading https://files.pythonhosted.org/packages/84/44/21a7fdd50841aaaef224b943f7d10df87e476e181bb926ccf859bcb53d48/graphviz-0.8.3-py2.py3-none-any.whl
Requirement already satisfied: requests<2.19.0,>=2.18.4 in /opt/conda/lib/python3.6/site-packages (from mxnet)
Requirement already satisfied: chardet<3.1.0,>=3.0.2 in /opt/conda/lib/python3.6/site-packages (from requests<2.19.0,>=2.18.4->mxnet)
Requirement already satisfied: idna<2.7,>=2.5 in /opt/conda/lib/python3.6/site-packages (from requests<2.19.0,>=2.18.4->mxnet)
Requirement already satisfied: urllib3<1.23,>=1.21.1 in /opt/conda/lib/python3.6/site-packages (from requests<2.19.0,>=2.18.4->mxnet)
Requirement already satisfied: certifi>=2017.4.17 in /opt/conda/lib/python3.6/site-packages (from requests<2.19.0,>=2.18.4->mxnet)
Installing collected packages: graphviz, mxnet
Successfully installed graphviz-0.8.3 mxnet-1.1.0.post0
You are using pip version 9.0.1, however version 10.0.1 is available.
You should consider upgrading via the 'pip install --upgrade pip' command.

import mxnet as mx
from mxnet import gluon
from scipy.stats import skew
from scipy.stats.stats import pearsonr
from mxnet import ndarray as nd
from mxnet import autograd
from mxnet import gluon

X_train,X_test=import_data()

all_X = pd.concat((X_train.loc[:, 'MSSubClass':'SaleCondition'],
                      X_test.loc[:, 'MSSubClass':'SaleCondition']))

numnic_feats = all_X.dtypes[all_X.dtypes!='object'].index
all_X[numnic_feats] = all_X[numnic_feats].apply(lambda x: (x-x.mean())/x.std())

all_X = pd.get_dummies(all_X, dummy_na=True)

all_X = all_X.fillna(all_X.mean())

num_train = X_train.shape[0]
X_train = all_X[:num_train].as_matrix()
X_test = all_X[num_train:].as_matrix()
y_train = y.as_matrix()

from mxnet import ndarray as nd
from mxnet import autograd
from mxnet import gluon

X_train = nd.array(X_train)
y_train = nd.array(y_train)
y_train.reshape((num_train, 1))

X_test = nd.array(X_test)

square_loss = gluon.loss.L2Loss()


def get_rmse_log(net, X_train,y_train):
    num_train = X_train.shape[0]
    clipped_preds = nd.clip(net(X_train), 1, float('inf'))
    return nd.sqrt(2*nd.sum(square_loss(nd.log(clipped_preds),nd.log(y_train))/num_train)).asscalar()

def get_net():
    net = gluon.nn.Sequential()
    with net.name_scope():
        net.add(gluon.nn.Dense(128, activation='relu'))
        net.add(gluon.nn.Dropout(0.01))
#         net.add(gluon.nn.Dense(32, activation='relu'))
#         net.add(gluon.nn.Dropout(0.2))
        net.add(gluon.nn.Dense(1))
    net.initialize()
    return net

%matplotlib inline
import matplotlib as mpl
mpl.rcParams['figure.dpi']= 120
import matplotlib.pyplot as plt

def train(net, X_train, y_train, X_test, y_test, epochs,
          verbose_epoch, learning_rate, weight_decay):
    train_loss = []
    if X_test is not None:
        test_loss = []
    batch_size = 100
    dataset_train = gluon.data.ArrayDataset(X_train, y_train)
    data_iter_train = gluon.data.DataLoader(
        dataset_train, batch_size,shuffle=True)
    trainer = gluon.Trainer(net.collect_params(), 'adam',
                            {'learning_rate': learning_rate,
                             'wd': weight_decay})
    net.collect_params().initialize(force_reinit=True)
    for epoch in range(epochs):
        for data, label in data_iter_train:
            with autograd.record():
                output = net(data)
                loss = square_loss(output, label)
            loss.backward()
            trainer.step(batch_size)

            cur_train_loss = get_rmse_log(net, X_train, y_train)
        if epoch > verbose_epoch:
            print("Epoch %d, train loss: %f" % (epoch, cur_train_loss))
        train_loss.append(cur_train_loss)
        if X_test is not None:
            cur_test_loss = get_rmse_log(net, X_test, y_test)
            test_loss.append(cur_test_loss)
    plt.plot(train_loss)
    plt.legend(['train'])
    if X_test is not None:
        plt.plot(test_loss)
        plt.legend(['train','test'])
    plt.show()
    if X_test is not None:
        return cur_train_loss, cur_test_loss
    else:
        return cur_train_loss

def k_fold_cross_valid(k, epochs, verbose_epoch, X_train, y_train,
                       learning_rate, weight_decay):
    assert k > 1
    fold_size = X_train.shape[0] // k
    train_loss_sum = 0.0
    test_loss_sum = 0.0
    for test_i in range(k):
        X_val_test = X_train[test_i * fold_size: (test_i + 1) * fold_size, :]
        y_val_test = y_train[test_i * fold_size: (test_i + 1) * fold_size]

        val_train_defined = False
        for i in range(k):
            if i != test_i:
                X_cur_fold = X_train[i * fold_size: (i + 1) * fold_size, :]
                y_cur_fold = y_train[i * fold_size: (i + 1) * fold_size]
                if not val_train_defined:
                    X_val_train = X_cur_fold
                    y_val_train = y_cur_fold
                    val_train_defined = True
                else:
                    X_val_train = nd.concat(X_val_train, X_cur_fold, dim=0)
                    y_val_train = nd.concat(y_val_train, y_cur_fold, dim=0)
        net = get_net()
        train_loss, test_loss = train(
            net, X_val_train, y_val_train, X_val_test, y_val_test,epochs, verbose_epoch, learning_rate, weight_decay)
        train_loss_sum += train_loss
        print("Test loss: %f" % test_loss)
        test_loss_sum += test_loss
    return train_loss_sum / k, test_loss_sum / k

k = 50
epochs = 100#80#100#80#100
verbose_epoch = 35#75#95
learning_rate = 0.01 #0.001 0.24
weight_decay = 130 #200#80 #100.0 #20
train_loss, test_loss = k_fold_cross_valid(k, epochs, verbose_epoch, X_train,
                                           y_train, learning_rate, weight_decay)
print("%d-fold validation: Avg train loss: %f, Avg test loss: %f" %
      (k, train_loss, test_loss))

Epoch 36, train loss: 0.358370
Epoch 37, train loss: 0.329547
Epoch 38, train loss: 0.304156
Epoch 39, train loss: 0.281300
Epoch 40, train loss: 0.260853
Epoch 41, train loss: 0.243371
Epoch 42, train loss: 0.228948
Epoch 43, train loss: 0.216911
Epoch 44, train loss: 0.206755
Epoch 45, train loss: 0.198232
Epoch 46, train loss: 0.191778
Epoch 47, train loss: 0.186597
Epoch 48, train loss: 0.182579
Epoch 49, train loss: 0.179636
Epoch 50, train loss: 0.177325
Epoch 51, train loss: 0.175631
Epoch 52, train loss: 0.174307
Epoch 53, train loss: 0.173259
Epoch 54, train loss: 0.172454
Epoch 55, train loss: 0.171768
Epoch 56, train loss: 0.171191
Epoch 57, train loss: 0.170652
Epoch 58, train loss: 0.170196
Epoch 59, train loss: 0.169767
Epoch 60, train loss: 0.169482
Epoch 61, train loss: 0.169156
Epoch 62, train loss: 0.168901
Epoch 63, train loss: 0.168714
Epoch 64, train loss: 0.168626
Epoch 65, train loss: 0.168415
Epoch 66, train loss: 0.168280
Epoch 67, train loss: 0.168266
Epoch 68, train loss: 0.168301
Epoch 69, train loss: 0.168384
Epoch 70, train loss: 0.168341
Epoch 71, train loss: 0.168410
Epoch 72, train loss: 0.168418
Epoch 73, train loss: 0.168599
Epoch 74, train loss: 0.168507
Epoch 75, train loss: 0.168739
Epoch 76, train loss: 0.168638
Epoch 77, train loss: 0.168682
Epoch 78, train loss: 0.168664
Epoch 79, train loss: 0.168227
Epoch 80, train loss: 0.168469
Epoch 81, train loss: 0.168285
Epoch 82, train loss: 0.168195
Epoch 83, train loss: 0.168076
Epoch 84, train loss: 0.167744
Epoch 85, train loss: 0.167307
Epoch 86, train loss: 0.167071
Epoch 87, train loss: 0.166384
Epoch 88, train loss: 0.166117
Epoch 89, train loss: 0.165722
Epoch 90, train loss: 0.165048
Epoch 91, train loss: 0.164526
Epoch 92, train loss: 0.163991
Epoch 93, train loss: 0.163025
Epoch 94, train loss: 0.162141
Epoch 95, train loss: 0.161371
Epoch 96, train loss: 0.160612
Epoch 97, train loss: 0.159548
Epoch 98, train loss: 0.158722
Epoch 99, train loss: 0.157808

Test loss: 0.140869
Epoch 36, train loss: 0.264208
Epoch 37, train loss: 0.244609
Epoch 38, train loss: 0.227945
Epoch 39, train loss: 0.214140
Epoch 40, train loss: 0.202930
Epoch 41, train loss: 0.194328
Epoch 42, train loss: 0.187346
Epoch 43, train loss: 0.182347
Epoch 44, train loss: 0.178700
Epoch 45, train loss: 0.176210
Epoch 46, train loss: 0.174308
Epoch 47, train loss: 0.172905
Epoch 48, train loss: 0.171876
Epoch 49, train loss: 0.170997
Epoch 50, train loss: 0.170306
Epoch 51, train loss: 0.169683
Epoch 52, train loss: 0.169178
Epoch 53, train loss: 0.168726
Epoch 54, train loss: 0.168284
Epoch 55, train loss: 0.167926
Epoch 56, train loss: 0.167705
Epoch 57, train loss: 0.167476
Epoch 58, train loss: 0.167378
Epoch 59, train loss: 0.167134
Epoch 60, train loss: 0.167086
Epoch 61, train loss: 0.167210
Epoch 62, train loss: 0.167235
Epoch 63, train loss: 0.167062
Epoch 64, train loss: 0.167313
Epoch 65, train loss: 0.167584
Epoch 66, train loss: 0.167431
Epoch 67, train loss: 0.167629
Epoch 68, train loss: 0.167705
Epoch 69, train loss: 0.167675
Epoch 70, train loss: 0.167822
Epoch 71, train loss: 0.167937
Epoch 72, train loss: 0.167755
Epoch 73, train loss: 0.167722
Epoch 74, train loss: 0.167583
Epoch 75, train loss: 0.167337
Epoch 76, train loss: 0.166973
Epoch 77, train loss: 0.166700
Epoch 78, train loss: 0.166719
Epoch 79, train loss: 0.166273
Epoch 80, train loss: 0.165787
Epoch 81, train loss: 0.165396
Epoch 82, train loss: 0.164940
Epoch 83, train loss: 0.164476
Epoch 84, train loss: 0.163614
Epoch 85, train loss: 0.162581
Epoch 86, train loss: 0.161665
Epoch 87, train loss: 0.160658
Epoch 88, train loss: 0.160159
Epoch 89, train loss: 0.158706
Epoch 90, train loss: 0.157462
Epoch 91, train loss: 0.156222
Epoch 92, train loss: 0.154904
Epoch 93, train loss: 0.153871
Epoch 94, train loss: 0.152790
Epoch 95, train loss: 0.151638
Epoch 96, train loss: 0.150837
Epoch 97, train loss: 0.149618
Epoch 98, train loss: 0.148630
Epoch 99, train loss: 0.147771

Test loss: 0.190945
Epoch 36, train loss: 0.267124
Epoch 37, train loss: 0.246855
Epoch 38, train loss: 0.229879
Epoch 39, train loss: 0.216016
Epoch 40, train loss: 0.204692
Epoch 41, train loss: 0.195842
Epoch 42, train loss: 0.189142
Epoch 43, train loss: 0.184189
Epoch 44, train loss: 0.180288
Epoch 45, train loss: 0.177445
Epoch 46, train loss: 0.175363
Epoch 47, train loss: 0.173874
Epoch 48, train loss: 0.172679
Epoch 49, train loss: 0.171772
Epoch 50, train loss: 0.171068
Epoch 51, train loss: 0.170505
Epoch 52, train loss: 0.170004
Epoch 53, train loss: 0.169629
Epoch 54, train loss: 0.169343
Epoch 55, train loss: 0.169045
Epoch 56, train loss: 0.168902
Epoch 57, train loss: 0.168848
Epoch 58, train loss: 0.168937
Epoch 59, train loss: 0.168878
Epoch 60, train loss: 0.169028
Epoch 61, train loss: 0.169277
Epoch 62, train loss: 0.169507
Epoch 63, train loss: 0.169763
Epoch 64, train loss: 0.170024
Epoch 65, train loss: 0.170285
Epoch 66, train loss: 0.170609
Epoch 67, train loss: 0.171271
Epoch 68, train loss: 0.171523
Epoch 69, train loss: 0.171807
Epoch 70, train loss: 0.172101
Epoch 71, train loss: 0.172324
Epoch 72, train loss: 0.172671
Epoch 73, train loss: 0.172784
Epoch 74, train loss: 0.173078
Epoch 75, train loss: 0.173085
Epoch 76, train loss: 0.173282
Epoch 77, train loss: 0.173501
Epoch 78, train loss: 0.173401
Epoch 79, train loss: 0.173387
Epoch 80, train loss: 0.173016
Epoch 81, train loss: 0.172710
Epoch 82, train loss: 0.172166
Epoch 83, train loss: 0.172064
Epoch 84, train loss: 0.171664
Epoch 85, train loss: 0.170945
Epoch 86, train loss: 0.170266
Epoch 87, train loss: 0.169716
Epoch 88, train loss: 0.169120
Epoch 89, train loss: 0.168491
Epoch 90, train loss: 0.167719
Epoch 91, train loss: 0.167067
Epoch 92, train loss: 0.166179
Epoch 93, train loss: 0.165225
Epoch 94, train loss: 0.164257
Epoch 95, train loss: 0.163442
Epoch 96, train loss: 0.162457
Epoch 97, train loss: 0.161400
Epoch 98, train loss: 0.160597
Epoch 99, train loss: 0.159378

Test loss: 0.190603
Epoch 36, train loss: 0.273645
Epoch 37, train loss: 0.252416
Epoch 38, train loss: 0.234936
Epoch 39, train loss: 0.220489
Epoch 40, train loss: 0.208710
Epoch 41, train loss: 0.199074
Epoch 42, train loss: 0.192034
Epoch 43, train loss: 0.186522
Epoch 44, train loss: 0.182302
Epoch 45, train loss: 0.179208
Epoch 46, train loss: 0.176844
Epoch 47, train loss: 0.175163
Epoch 48, train loss: 0.173855
Epoch 49, train loss: 0.172844
Epoch 50, train loss: 0.172074
Epoch 51, train loss: 0.171444
Epoch 52, train loss: 0.170913
Epoch 53, train loss: 0.170477
Epoch 54, train loss: 0.170108
Epoch 55, train loss: 0.169785
Epoch 56, train loss: 0.169509
Epoch 57, train loss: 0.169336
Epoch 58, train loss: 0.169193
Epoch 59, train loss: 0.169254
Epoch 60, train loss: 0.169178
Epoch 61, train loss: 0.169215
Epoch 62, train loss: 0.169297
Epoch 63, train loss: 0.169325
Epoch 64, train loss: 0.169711
Epoch 65, train loss: 0.169907
Epoch 66, train loss: 0.170293
Epoch 67, train loss: 0.170751
Epoch 68, train loss: 0.170831
Epoch 69, train loss: 0.171007
Epoch 70, train loss: 0.171588
Epoch 71, train loss: 0.171934
Epoch 72, train loss: 0.171955
Epoch 73, train loss: 0.172107
Epoch 74, train loss: 0.172334
Epoch 75, train loss: 0.172215
Epoch 76, train loss: 0.172334
Epoch 77, train loss: 0.171943
Epoch 78, train loss: 0.171994
Epoch 79, train loss: 0.171767
Epoch 80, train loss: 0.171410
Epoch 81, train loss: 0.170318
Epoch 82, train loss: 0.169657
Epoch 83, train loss: 0.169487
Epoch 84, train loss: 0.168771
Epoch 85, train loss: 0.168205
Epoch 86, train loss: 0.167218
Epoch 87, train loss: 0.166725
Epoch 88, train loss: 0.166045
Epoch 89, train loss: 0.164897
Epoch 90, train loss: 0.164499
Epoch 91, train loss: 0.163571
Epoch 92, train loss: 0.162622
Epoch 93, train loss: 0.161339
Epoch 94, train loss: 0.160432
Epoch 95, train loss: 0.159999
Epoch 96, train loss: 0.159022
Epoch 97, train loss: 0.158091
Epoch 98, train loss: 0.157329
Epoch 99, train loss: 0.156679

Test loss: 0.173561
Epoch 36, train loss: 0.265279
Epoch 37, train loss: 0.245101
Epoch 38, train loss: 0.228695
Epoch 39, train loss: 0.215127
Epoch 40, train loss: 0.204086
Epoch 41, train loss: 0.195763
Epoch 42, train loss: 0.189128
Epoch 43, train loss: 0.184137
Epoch 44, train loss: 0.180524
Epoch 45, train loss: 0.177690
Epoch 46, train loss: 0.175740
Epoch 47, train loss: 0.174208
Epoch 48, train loss: 0.173083
Epoch 49, train loss: 0.172179
Epoch 50, train loss: 0.171495
Epoch 51, train loss: 0.170944
Epoch 52, train loss: 0.170377
Epoch 53, train loss: 0.170011
Epoch 54, train loss: 0.169665
Epoch 55, train loss: 0.169378
Epoch 56, train loss: 0.169133
Epoch 57, train loss: 0.169036
Epoch 58, train loss: 0.168782
Epoch 59, train loss: 0.168709
Epoch 60, train loss: 0.168785
Epoch 61, train loss: 0.168810
Epoch 62, train loss: 0.168805
Epoch 63, train loss: 0.168835
Epoch 64, train loss: 0.168878
Epoch 65, train loss: 0.169193
Epoch 66, train loss: 0.169227
Epoch 67, train loss: 0.169245
Epoch 68, train loss: 0.169561
Epoch 69, train loss: 0.169529
Epoch 70, train loss: 0.169516
Epoch 71, train loss: 0.169371
Epoch 72, train loss: 0.169330
Epoch 73, train loss: 0.169298
Epoch 74, train loss: 0.169057
Epoch 75, train loss: 0.168943
Epoch 76, train loss: 0.168509
Epoch 77, train loss: 0.168267
Epoch 78, train loss: 0.167600
Epoch 79, train loss: 0.167156
Epoch 80, train loss: 0.166645
Epoch 81, train loss: 0.165861
Epoch 82, train loss: 0.165372
Epoch 83, train loss: 0.164906
Epoch 84, train loss: 0.164218
Epoch 85, train loss: 0.163475
Epoch 86, train loss: 0.162663
Epoch 87, train loss: 0.162144
Epoch 88, train loss: 0.161552
Epoch 89, train loss: 0.160771
Epoch 90, train loss: 0.159926
Epoch 91, train loss: 0.159248
Epoch 92, train loss: 0.158681
Epoch 93, train loss: 0.157823
Epoch 94, train loss: 0.156875
Epoch 95, train loss: 0.156089
Epoch 96, train loss: 0.155575
Epoch 97, train loss: 0.154678
Epoch 98, train loss: 0.153843
Epoch 99, train loss: 0.153167

Test loss: 0.148445
Epoch 36, train loss: 0.254571
Epoch 37, train loss: 0.236752
Epoch 38, train loss: 0.221799
Epoch 39, train loss: 0.209882
Epoch 40, train loss: 0.200598
Epoch 41, train loss: 0.193127
Epoch 42, train loss: 0.187378
Epoch 43, train loss: 0.183208
Epoch 44, train loss: 0.180053
Epoch 45, train loss: 0.177665
Epoch 46, train loss: 0.175897
Epoch 47, train loss: 0.174478
Epoch 48, train loss: 0.173434
Epoch 49, train loss: 0.172600
Epoch 50, train loss: 0.171897
Epoch 51, train loss: 0.171342
Epoch 52, train loss: 0.170879
Epoch 53, train loss: 0.170467
Epoch 54, train loss: 0.170152
Epoch 55, train loss: 0.169941
Epoch 56, train loss: 0.169754
Epoch 57, train loss: 0.169635
Epoch 58, train loss: 0.169607
Epoch 59, train loss: 0.169664
Epoch 60, train loss: 0.169815
Epoch 61, train loss: 0.169781
Epoch 62, train loss: 0.169928
Epoch 63, train loss: 0.170118
Epoch 64, train loss: 0.170463
Epoch 65, train loss: 0.170681
Epoch 66, train loss: 0.171009
Epoch 67, train loss: 0.171263
Epoch 68, train loss: 0.171901
Epoch 69, train loss: 0.172124
Epoch 70, train loss: 0.172274
Epoch 71, train loss: 0.172448
Epoch 72, train loss: 0.172890
Epoch 73, train loss: 0.172989
Epoch 74, train loss: 0.173297
Epoch 75, train loss: 0.173366
Epoch 76, train loss: 0.173207
Epoch 77, train loss: 0.173098
Epoch 78, train loss: 0.172752
Epoch 79, train loss: 0.172677
Epoch 80, train loss: 0.172595
Epoch 81, train loss: 0.172219
Epoch 82, train loss: 0.171680
Epoch 83, train loss: 0.171281
Epoch 84, train loss: 0.170442
Epoch 85, train loss: 0.169488
Epoch 86, train loss: 0.169250
Epoch 87, train loss: 0.168372
Epoch 88, train loss: 0.167834
Epoch 89, train loss: 0.166859
Epoch 90, train loss: 0.165998
Epoch 91, train loss: 0.164927
Epoch 92, train loss: 0.164132
Epoch 93, train loss: 0.163411
Epoch 94, train loss: 0.162900
Epoch 95, train loss: 0.161870
Epoch 96, train loss: 0.161251
Epoch 97, train loss: 0.160105
Epoch 98, train loss: 0.159418
Epoch 99, train loss: 0.158560

Test loss: 0.126427
Epoch 36, train loss: 0.292918
Epoch 37, train loss: 0.270132
Epoch 38, train loss: 0.250267
Epoch 39, train loss: 0.233007
Epoch 40, train loss: 0.219400
Epoch 41, train loss: 0.208049
Epoch 42, train loss: 0.198942
Epoch 43, train loss: 0.191799
Epoch 44, train loss: 0.186212
Epoch 45, train loss: 0.182173
Epoch 46, train loss: 0.179162
Epoch 47, train loss: 0.176740
Epoch 48, train loss: 0.175066
Epoch 49, train loss: 0.173699
Epoch 50, train loss: 0.172633
Epoch 51, train loss: 0.171805
Epoch 52, train loss: 0.171105
Epoch 53, train loss: 0.170566
Epoch 54, train loss: 0.170030
Epoch 55, train loss: 0.169626
Epoch 56, train loss: 0.169295
Epoch 57, train loss: 0.168898
Epoch 58, train loss: 0.168711
Epoch 59, train loss: 0.168451
Epoch 60, train loss: 0.168304
Epoch 61, train loss: 0.168354
Epoch 62, train loss: 0.168192
Epoch 63, train loss: 0.168259
Epoch 64, train loss: 0.168510
Epoch 65, train loss: 0.168605
Epoch 66, train loss: 0.168641
Epoch 67, train loss: 0.168810
Epoch 68, train loss: 0.168970
Epoch 69, train loss: 0.169289
Epoch 70, train loss: 0.169463
Epoch 71, train loss: 0.169833
Epoch 72, train loss: 0.170005
Epoch 73, train loss: 0.170122
Epoch 74, train loss: 0.170213
Epoch 75, train loss: 0.170155
Epoch 76, train loss: 0.170224
Epoch 77, train loss: 0.170386
Epoch 78, train loss: 0.170118
Epoch 79, train loss: 0.170216
Epoch 80, train loss: 0.169969
Epoch 81, train loss: 0.169738
Epoch 82, train loss: 0.169438
Epoch 83, train loss: 0.168955
Epoch 84, train loss: 0.168751
Epoch 85, train loss: 0.168574
Epoch 86, train loss: 0.167916
Epoch 87, train loss: 0.167562
Epoch 88, train loss: 0.166956
Epoch 89, train loss: 0.166104
Epoch 90, train loss: 0.165800
Epoch 91, train loss: 0.165086
Epoch 92, train loss: 0.164388
Epoch 93, train loss: 0.163435
Epoch 94, train loss: 0.163043
Epoch 95, train loss: 0.162289
Epoch 96, train loss: 0.161542
Epoch 97, train loss: 0.160862
Epoch 98, train loss: 0.159995
Epoch 99, train loss: 0.159472

Test loss: 0.159065
Epoch 36, train loss: 0.292998
Epoch 37, train loss: 0.268990
Epoch 38, train loss: 0.248293
Epoch 39, train loss: 0.231621
Epoch 40, train loss: 0.217530
Epoch 41, train loss: 0.206223
Epoch 42, train loss: 0.197537
Epoch 43, train loss: 0.190819
Epoch 44, train loss: 0.185533
Epoch 45, train loss: 0.181788
Epoch 46, train loss: 0.178885
Epoch 47, train loss: 0.176664
Epoch 48, train loss: 0.175000
Epoch 49, train loss: 0.173708
Epoch 50, train loss: 0.172631
Epoch 51, train loss: 0.171786
Epoch 52, train loss: 0.171071
Epoch 53, train loss: 0.170434
Epoch 54, train loss: 0.169877
Epoch 55, train loss: 0.169386
Epoch 56, train loss: 0.168949
Epoch 57, train loss: 0.168565
Epoch 58, train loss: 0.168150
Epoch 59, train loss: 0.167883
Epoch 60, train loss: 0.167657
Epoch 61, train loss: 0.167299
Epoch 62, train loss: 0.167167
Epoch 63, train loss: 0.167119
Epoch 64, train loss: 0.166881
Epoch 65, train loss: 0.166678
Epoch 66, train loss: 0.166592
Epoch 67, train loss: 0.166494
Epoch 68, train loss: 0.166222
Epoch 69, train loss: 0.166096
Epoch 70, train loss: 0.165934
Epoch 71, train loss: 0.165800
Epoch 72, train loss: 0.165329
Epoch 73, train loss: 0.165014
Epoch 74, train loss: 0.164549
Epoch 75, train loss: 0.164194
Epoch 76, train loss: 0.163766
Epoch 77, train loss: 0.163291
Epoch 78, train loss: 0.162880
Epoch 79, train loss: 0.162220
Epoch 80, train loss: 0.161684
Epoch 81, train loss: 0.161166
Epoch 82, train loss: 0.160534
Epoch 83, train loss: 0.159822
Epoch 84, train loss: 0.159135
Epoch 85, train loss: 0.158559
Epoch 86, train loss: 0.158080
Epoch 87, train loss: 0.157206
Epoch 88, train loss: 0.156699
Epoch 89, train loss: 0.155876
Epoch 90, train loss: 0.155251
Epoch 91, train loss: 0.154578
Epoch 92, train loss: 0.153733
Epoch 93, train loss: 0.152925
Epoch 94, train loss: 0.152265
Epoch 95, train loss: 0.151495
Epoch 96, train loss: 0.150724
Epoch 97, train loss: 0.150214
Epoch 98, train loss: 0.149322
Epoch 99, train loss: 0.148771

Test loss: 0.165271
Epoch 36, train loss: 0.250949
Epoch 37, train loss: 0.233455
Epoch 38, train loss: 0.219188
Epoch 39, train loss: 0.207545
Epoch 40, train loss: 0.198473
Epoch 41, train loss: 0.191321
Epoch 42, train loss: 0.186009
Epoch 43, train loss: 0.182190
Epoch 44, train loss: 0.179151
Epoch 45, train loss: 0.177011
Epoch 46, train loss: 0.175339
Epoch 47, train loss: 0.174100
Epoch 48, train loss: 0.173160
Epoch 49, train loss: 0.172379
Epoch 50, train loss: 0.171742
Epoch 51, train loss: 0.171249
Epoch 52, train loss: 0.170781
Epoch 53, train loss: 0.170417
Epoch 54, train loss: 0.170094
Epoch 55, train loss: 0.169793
Epoch 56, train loss: 0.169663
Epoch 57, train loss: 0.169627
Epoch 58, train loss: 0.169583
Epoch 59, train loss: 0.169443
Epoch 60, train loss: 0.169544
Epoch 61, train loss: 0.169571
Epoch 62, train loss: 0.169689
Epoch 63, train loss: 0.169907
Epoch 64, train loss: 0.170006
Epoch 65, train loss: 0.170192
Epoch 66, train loss: 0.170578
Epoch 67, train loss: 0.170896
Epoch 68, train loss: 0.171122
Epoch 69, train loss: 0.171017
Epoch 70, train loss: 0.171141
Epoch 71, train loss: 0.171435
Epoch 72, train loss: 0.171385
Epoch 73, train loss: 0.171184
Epoch 74, train loss: 0.171251
Epoch 75, train loss: 0.171271
Epoch 76, train loss: 0.171072
Epoch 77, train loss: 0.170670
Epoch 78, train loss: 0.170269
Epoch 79, train loss: 0.170145
Epoch 80, train loss: 0.169753
Epoch 81, train loss: 0.169300
Epoch 82, train loss: 0.169110
Epoch 83, train loss: 0.168266
Epoch 84, train loss: 0.167807
Epoch 85, train loss: 0.167228
Epoch 86, train loss: 0.166356
Epoch 87, train loss: 0.165529
Epoch 88, train loss: 0.164774
Epoch 89, train loss: 0.164271
Epoch 90, train loss: 0.163834
Epoch 91, train loss: 0.162610
Epoch 92, train loss: 0.162001
Epoch 93, train loss: 0.161208
Epoch 94, train loss: 0.160338
Epoch 95, train loss: 0.159803
Epoch 96, train loss: 0.159031
Epoch 97, train loss: 0.158370
Epoch 98, train loss: 0.157515
Epoch 99, train loss: 0.157273

Test loss: 0.118897
Epoch 36, train loss: 0.302814
Epoch 37, train loss: 0.278594
Epoch 38, train loss: 0.257886
Epoch 39, train loss: 0.240270
Epoch 40, train loss: 0.225310
Epoch 41, train loss: 0.212968
Epoch 42, train loss: 0.202917
Epoch 43, train loss: 0.195051
Epoch 44, train loss: 0.188898
Epoch 45, train loss: 0.183969
Epoch 46, train loss: 0.180402
Epoch 47, train loss: 0.177718
Epoch 48, train loss: 0.175778
Epoch 49, train loss: 0.174239
Epoch 50, train loss: 0.173081
Epoch 51, train loss: 0.172178
Epoch 52, train loss: 0.171424
Epoch 53, train loss: 0.170790
Epoch 54, train loss: 0.170317
Epoch 55, train loss: 0.169876
Epoch 56, train loss: 0.169422
Epoch 57, train loss: 0.169102
Epoch 58, train loss: 0.168777
Epoch 59, train loss: 0.168538
Epoch 60, train loss: 0.168261
Epoch 61, train loss: 0.168162
Epoch 62, train loss: 0.168065
Epoch 63, train loss: 0.167951
Epoch 64, train loss: 0.167922
Epoch 65, train loss: 0.167754
Epoch 66, train loss: 0.167708
Epoch 67, train loss: 0.167814
Epoch 68, train loss: 0.167984
Epoch 69, train loss: 0.168065
Epoch 70, train loss: 0.168033
Epoch 71, train loss: 0.168039
Epoch 72, train loss: 0.168029
Epoch 73, train loss: 0.168147
Epoch 74, train loss: 0.168086
Epoch 75, train loss: 0.168157
Epoch 76, train loss: 0.167955
Epoch 77, train loss: 0.167929
Epoch 78, train loss: 0.167759
Epoch 79, train loss: 0.167246
Epoch 80, train loss: 0.167123
Epoch 81, train loss: 0.166591
Epoch 82, train loss: 0.166147
Epoch 83, train loss: 0.165975
Epoch 84, train loss: 0.165286
Epoch 85, train loss: 0.164982
Epoch 86, train loss: 0.164045
Epoch 87, train loss: 0.163376
Epoch 88, train loss: 0.162535
Epoch 89, train loss: 0.162105
Epoch 90, train loss: 0.161528
Epoch 91, train loss: 0.160699
Epoch 92, train loss: 0.159839
Epoch 93, train loss: 0.159214
Epoch 94, train loss: 0.158296
Epoch 95, train loss: 0.157426
Epoch 96, train loss: 0.156954
Epoch 97, train loss: 0.155929
Epoch 98, train loss: 0.155280
Epoch 99, train loss: 0.154294

By looking at the average values,our stacked model is the best.

Conclusions and Remarks

The best model seems to be the stacked one.
Other ideas would be to following a bayesian inference approach and implementing a gaussian process which we will not do for the lack of time and since we have a large dataset and inverting it takes time.
Deepnets could produce better results but need a lot of tweaking and experimenting.
Getting a lower RMSE needs experimenting with various things. We have no perfect data not perfect data that can give us perfect results,we just have to try.
Some cells were not run again since they take too much time.

References

A big thanks goes to Prof. Pietro Michiardi and the teaching assistant for answering my questions and guiding me.

The following is al ist of references we used:
https://machinelearningmastery.com/discover-feature-engineering-how-to-engineer-features-and-how-to-get-good-at-it/
https://machinelearningmastery.com/an-introduction-to-feature-selection/
http://jmlr.csail.mit.edu/papers/volume3/guyon03a/guyon03a.pdf
https://www.analyticsvidhya.com/blog/2016/01/complete-tutorial-ridge-lasso-regression-python/
https://codingstartups.com/practical-machine-learning-ridge-regression-vs-lasso/
https://machinelearningmastery.com/tutorial-first-neural-network-python-keras/https://codingstartups.com/practical-machine-learning-ridge-regression-vs-lasso/
https://www.kaggle.com/erikbruin/house-prices-lasso-xgboost-and-a-detailed-eda
https://www.kaggle.com/serigne/stacked-regressions-top-4-on-leaderboard
http://blog.kaggle.com/2017/06/15/stacking-made-easy-an-introduction-to-stacknet-by-competitions-grandmaster-marios-michailidis-kazanova/
https://www.kaggle.com/agehsbarg/top-10-0-10943-stacking-mice-and-brutal-force
https://www.kaggle.com/humananalog/xgboost-lasso
https://www.kaggle.com/apapiu/regularized-linear-models
http://ww2.amstat.org/publications/jse/v19n3/Decock/DataDocumentation.txt
https://www.kaggle.com/harunshimanto/house-price-prediction

Number of variables	81
Number of observations	1200
Total Missing (%)	5.9%
Total size in memory	759.5 KiB
Average record size in memory	648.1 B

Numeric	38
Categorical	43
Boolean	0
Date	0
Text (Unique)	0
Rejected	0
Unsupported	0

Distinct count	676
Unique (%)	56.3%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Mean	1157.4
Minimum	334
Maximum	3228
Zeros (%)	0.0%

Minimum	334
5-th percentile	672
Q1	882
Median	1087
Q3	1390.2
95-th percentile	1803.4
Maximum	3228
Range	2894
Interquartile range	508.25

Challenge: What is the future selling price of a home?

Approach #1

Descriptive Statistics of the data

Overview

Variables

Correlations

Sample

Data Cleaning and Feature Selection

Testing Different Models and Comparing

Remarks and Notes

Approach #2

Data Cleaning

Overview

Variables

Correlations

Sample

Feature Engineering and Selection

Training Time

Approach #3

Feature Engineering

Training

Approach #4

Conclusions and Remarks

References

	Id	MSSubClass	MSZoning	LotFrontage	LotArea	Street	Alley	LotShape	LandContour	Utilities	...	PoolQC	Fence	MiscFeature	MoSold	YrSold	SaleType	SaleCondition	SalePrice
0	1	60	RL	65.0	8450	Pave	NaN	Reg	Lvl	AllPub	...	NaN	NaN	NaN	2	2008	WD	Normal	208500
1	2	20	RL	80.0	9600	Pave	NaN	Reg	Lvl	AllPub	...	NaN	NaN	NaN	5	2007	WD	Normal	181500
2	3	60	RL	68.0	11250	Pave	NaN	IR1	Lvl	AllPub	...	NaN	NaN	NaN	9	2008	WD	Normal	223500
3	4	70	RL	60.0	9550	Pave	NaN	IR1	Lvl	AllPub	...	NaN	NaN	NaN	2	2006	WD	Abnorml	140000
4	5	60	RL	84.0	14260	Pave	NaN	IR1	Lvl	AllPub	...	NaN	NaN	NaN	12	2008	WD	Normal	250000

	Id	MSSubClass	MSZoning	LotFrontage	LotArea	Street	Alley	LotShape	LandContour	Utilities	...	PoolQC	Fence	MiscFeature	MoSold	YrSold	SaleType	SaleCondition
0	1201	20	RL	71.0	9353	Pave	NaN	Reg	Lvl	AllPub	...	NaN	NaN	Shed	7	2006	Oth	Abnorml
1	1202	60	RL	80.0	10400	Pave	NaN	Reg	Lvl	AllPub	...	NaN	NaN	NaN	3	2009	WD	Normal
2	1203	50	RM	50.0	6000	Pave	NaN	Reg	Lvl	AllPub	...	NaN	NaN	NaN	5	2009	WD	Normal
3	1204	20	RL	75.0	9750	Pave	NaN	Reg	Lvl	AllPub	...	NaN	NaN	NaN	10	2009	WD	Normal
4	1205	20	RL	78.0	10140	Pave	NaN	Reg	Lvl	AllPub	...	NaN	MnPrv	NaN	7	2006	WD	Normal

Standard deviation	375.24
Coef of variation	0.3242
Kurtosis	1.7089
Mean	1157.4
MAD	296.97
Skewness	0.96693
Sum	1388917
Variance	140800
Memory size	9.5 KiB

Value	Count	Frequency (%)
864	18	1.5%
1040	13	1.1%
912	11	0.9%
894	10	0.8%
848	9	0.8%
672	8	0.7%
630	7	0.6%
936	7	0.6%
816	7	0.6%
832	6	0.5%
Other values (666)	1104	92.0%

Distinct count	361
Unique (%)	30.1%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Standard deviation	437.04
Coef of variation	1.2668
Kurtosis	-0.48487
Mean	344.99
MAD	396.1
Skewness	0.83545
Sum	413992
Variance	191000
Memory size	9.5 KiB

Standard deviation	29.991
Coef of variation	8.2092
Kurtosis	122.08
Mean	3.6533
MAD	7.1727
Skewness	10.122
Sum	4384
Variance	899.47
Memory size	9.5 KiB

Distinct count	3
Unique (%)	0.2%
Missing (%)	93.8%
Missing (n)	1125

Standard deviation	0.8192
Coef of variation	0.28668
Kurtosis	2.3648
Mean	2.8575
MAD	0.58259
Skewness	0.27726
Sum	3429
Variance	0.67109
Memory size	9.5 KiB

Distinct count	567
Unique (%)	47.2%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Value	Count	Frequency (%)
334	1	0.1%
372	1	0.1%
438	1	0.1%
480	1	0.1%
483	6	0.5%

Value	Count	Frequency (%)
2515	1	0.1%
2524	1	0.1%
2898	1	0.1%
3138	1	0.1%
3228	1	0.1%

Value	Count	Frequency (%)
0	685	57.1%
728	8	0.7%
504	8	0.7%
720	7	0.6%
546	7	0.6%
672	6	0.5%
689	5	0.4%
840	5	0.4%
780	5	0.4%
756	5	0.4%
Other values (351)	459	38.2%

Value	Count	Frequency (%)
0	685	57.1%
110	1	0.1%
167	1	0.1%
213	1	0.1%
220	1	0.1%

Standard deviation	439.99
Coef of variation	0.98899
Kurtosis	0.075457
Mean	444.89
MAD	366.63
Skewness	0.80561
Sum	533864
Variance	193590
Memory size	9.5 KiB

Distinct count	116
Unique (%)	9.7%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Standard deviation	158.93
Coef of variation	3.5115
Kurtosis	21.374
Mean	45.26
MAD	80.503
Skewness	4.35
Sum	54312
Variance	25259
Memory size	9.5 KiB

Distinct count	4
Unique (%)	0.3%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Standard deviation	0.52034
Coef of variation	1.234
Kurtosis	-0.71336
Mean	0.42167
MAD	0.49827
Skewness	0.63666
Sum	506
Variance	0.27076
Memory size	9.5 KiB

Standard deviation	0.24748
Coef of variation	4.0133
Kurtosis	15.43
Mean	0.061667
MAD	0.11593
Skewness	3.9756
Sum	74
Variance	0.061248
Memory size	9.5 KiB

Distinct count	691
Unique (%)	57.6%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Standard deviation	440.39
Coef of variation	0.78007
Kurtosis	0.58593
Mean	564.55
MAD	350.12
Skewness	0.95092
Sum	677464
Variance	193940
Memory size	9.5 KiB

Distinct count	104
Unique (%)	8.7%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Standard deviation	61.507
Coef of variation	2.7733
Kurtosis	10.979
Mean	22.178
MAD	37.93
Skewness	3.1286
Sum	26614
Variance	3783.2
Memory size	9.5 KiB

Standard deviation	0.64211
Coef of variation	1.0455
Kurtosis	-0.22256
Mean	0.61417
MAD	0.57732
Skewness	0.63788
Sum	737
Variance	0.41231
Memory size	9.5 KiB

Standard deviation	0.55217
Coef of variation	0.35376
Kurtosis	-0.85784
Mean	1.5608
MAD	0.52389
Skewness	0.070791
Sum	1873
Variance	0.30489
Memory size	9.5 KiB

Distinct count	402
Unique (%)	33.5%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0

Minimum	0
5-th percentile	0
Q1	334.5
Median	478
Q3	576
95-th percentile	852.05
Maximum	1390
Range	1390
Interquartile range	241.5

Standard deviation	212.72
Coef of variation	0.45011
Kurtosis	0.82872
Mean	472.6
MAD	159.76
Skewness	0.14444
Sum	567125
Variance	45251
Memory size	9.5 KiB

Standard deviation	0.74349
Coef of variation	0.42164
Kurtosis	0.20649
Mean	1.7633
MAD	0.58239
Skewness	-0.36896
Sum	2116
Variance	0.55278
Memory size	9.5 KiB

Distinct count	96
Unique (%)	8.0%
Missing (%)	5.6%
Missing (n)	67
Infinite (%)	0.0%
Infinite (n)	0

Minimum	1900
5-th percentile	1930
Q1	1961
Median	1980
Q3	2002
95-th percentile	2007
Maximum	2010
Range	110
Interquartile range	41

Standard deviation	24.813
Coef of variation	0.012542
Kurtosis	-0.41703
Mean	1978.4
MAD	20.99
Skewness	-0.6562
Sum	2241500
Variance	615.68
Memory size	9.5 KiB

Distinct count	761
Unique (%)	63.4%
Missing (%)	0.0%
Missing (n)	0
Infinite (%)	0.0%
Infinite (n)	0